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EQUATIONS OF MOTION OF SLUNG-LOAD SYSTEMS, INCLUDING 

MULTILIFT SYSTEMS 

Luigi S. Cicolani and Gerd Kanning 
Ames Research Center 


SUMMARY 

General simulation equations are derived for the 
rigid body motion of slung-load systems. This work is 
motivated by an interest in trajectory control for slung 
loads carried by two or more helicopters. An approxi- 
mation of these systems consists of several rigid bodies 
connected by straight-line cables or links. The sus- 
pension can be assumed elastic or inelastic. Equations 
for the general system are obtained from the Newton- 
Euler rigid-body equations with the introduction of 
generalized velocity coordinates. Three forms are ob- 
tained: two generalize previous case-specific results for 
single-helicopter systems with elastic and inelastic sus- 
pensions, respectively, and the third is a new formu- 
lation for inelastic suspensions. The latter is derived 
from the elastic suspension equations by choosing the 
generalized coordinates so that motion induced by ca- 
ble stretching is separated from motion with invariant 
cable lengths, and by then nulling the stretching coor- 
dinates to get a relation for the suspension forces. The 
result is computationally more efficient than the con- 
ventional formulation, is readily integrated with the 
elastic suspension formulation, and is easily applied 
to the complex dual-lift and multilift systems. Re- 
sults are given for two- helicopter systems; three con- 
figurations are included and these can be integrated 
in a single simulation. Equations are also given for 
some single-helicopter systems, for comparison with 
the previous literature, and for a multilift system. 
Equations for degenerate-body approximations (point 
masses, rigid rods) are also formulated and results are 
given for dual-lift and multilift systems. Finally, lin- 
earized equations of motion are given for general slung- 
load systems are presented along with results for the 
two-helicopter system with a spreader bar. 


1. INTRODUCTION 
Background 

Various actual and proposed slung-load systems are 
illustrated in figure 1. Single-helicopter slung-load op- 
erations with the load suspended by cables in various 
ways from a single attachment point have been com- 
mon since the 1950s. Such operations were further 
developed and extensively used during the Vietnam 
war. Later research during the period 1965-1975 for 
the Heavy Lift Helicopter was focused on the stabi- 
lization of difficult loads, such as the standard 8- by 8- 
by 20-ft cargo container (MILVAN), by means of sus- 
pensions with multiple attachment points and various 
control devices. 

The use of two or more helicopters has been pe- 
riodically proposed since the early success of single- 
helicopter operations using suspensions consisting of 
cables and spreader bars (systems 8-12 in fig. 1) 
(refs. 1-3). Dual-lift suspensions have received lim- 
ited flight testing, have been used to carry payloads in 
a few isolated commercial operations, and have been 
advocated as an alternative either to developing new 
helicopters with greater payloads than those of cur- 
rent helicopters (ref. 4), or to obtaining larger-capacity 
helicopters than those locally available in a given sit- 
uation. A significant obstacle to further operational 
development is the complexity of system motion and 
its guidance and control along any typical maneuver- 
ing flight path (ref. 5). Until recently, progress beyond 
the initial investigation of hover dynamics and control 
in references 6 and 7 has been hampered by the lack 
of realistic and comprehensive equations of motion for 
use in theoretical and simulation studies. Tractabil- 
ity of the equations for analysis and programming and 
computational efficiency become critical factors for the 
multilift systems. Whereas the slung-load systems can 
be viewed simply as a few rigid bodies connected by 
cables, considerable complexity of the equations of mo- 
tion (EOMs) arises in applying the classic methods of 


n = Number of rigid bodies 
c = Number of constraints (inelastic suspension) 
d = 6n-c = number of degrees of freedom 
m = Number of suspension lines 



Figure 1. Examples of slung-load systems. 



the previous slung-load literature to the multilift sys- 
tems when the cables are modeled as inelastic. 

The approach taken in this paper and in our ini- 
tial report (ref. 8) is to develop a systematic analyti- 
cal formulation for general slung-load systems and an- 
alytical devices tailored to such systems, which are 
readily applied to specific cases and yield tractable, 
efficient equations. The devices include special coordi- 
nates that represent the suspension constraints in an 
inelastic system, and the systematic use of the natural 
vectors and matrices of rigid-body mechanics in the 
applications work. 

An alternative approach that would circumvent the 
need for further analysis would be to use one of the 
commercially available computer programs for the dy- 
namics of general multibody systems. A third poten- 
tial approach would be to apply the previous deriva- 
tion techniques to multilift systems by using symbolic 
digital computations to circumvent the excessive labor 
and unfavorable error probabilities of extended hand 
analysis and programming. However, these alternative 
approaches are designed to accommodate a very large 
class of dynamic systems, and the results do not pro- 
vide the insight and computational efficiency achieved 
here by restricting the class to slung-load systems. 

Equations of Motion for Slung-Load Systems 

The slung-load systems shown in figure 1 are viewed 
here as members of a class of systems consisting of rigid 
bodies connected by massless straight-line links which 
can be either elastic or inelastic, and w f hich support 
only forces along the link. These systems are charac- 
terized by the mass, inertia, and aerodynamics of the 
rigid bodies, and the suspension’s attachment point 
locations, unloaded link lengths, and link elastic pa- 
rameters. The limitations of these class properties in 
representing the slung-load systems are as follows. The 
rigid-body assumption excludes flexible loads and heli- 
copter elastic modes; cable mass and aerodynamics are 
neglected; and cable stretching is neglected in the case 
of inelastic cables. Despite these limitations, the class 
properties are expected to suffice for trajectory-control 
studies in which only low-frequency phenomena are of 
interest. 

Previous derivations of the EOMs for single- 
helicopter slung-load systems have assumed either in- 
elastic cables (refs. 9-12), or elastic cables (refs. 13-15), 
or considered both cable models (ref. 16). In most of 


these works, the results are specific to particular sus- 
pensions, but references 13 and 15 account for a gen- 
eral set of elastic suspensions in which all cables con- 
nect two rigid bodies. Early work on the EOMs for 
multilift systems assumed inelastic cables and yielded 
only limited results. In reference 17, a general formula- 
tion for systems of point masses connected by inelastic 
cables is given, and reference 18 contains equations for 
the same approximation of the three-body, dual-lift 
system 9 of figure 1. In references 6 and 7 linearized 
hover equations are given for the four-body, dual-lift 
system 10 of figure 1, assuming a point mass load. 

In the above-cited literature, the formulation of the 
EOMs depends principally on their use for simulation 
or control analysis. For simulations, the suspension is 
usually assumed elastic, and the equations are formu- 
lated from the Newton-Euler equations for the rigid- 
body velocity coordinates as, for example, in refer- 
ences 13-15. This method can be readily applied to 
most of the systems shown in figure 1, including the 
dual-lift and multilift systems. The forces and mo- 
ments applied to the configuration by the cables ap- 
pear explicitly in the equations, and are calculated 
from the system geometry and from a simple spring 
model of cable stretching. 

In practice, however, cables are relatively stiff, so the 
rigid-body motion with inelastic cables differs very lit- 
tle from that with elastic cables, and a reduced-order 
inelastic-cable model is of interest in trajectory-control 
analysis. Equations of motion for inelastic cables have 
usually been derived from Lagrange’s equations for 
general dynamical systems. Equations are given for a 
minimal set of suspension-specific, generalized position 
coordinates, and the nonworking suspension forces 
are eliminated, consistent with d’Alembert’s princi- 
ple. These equations require the inverse of a large 
system matrix and contain lengthy second-order ve- 
locity terms in exchange for the reduced system order 
and eliminated suspension forces. These features ob- 
struct the use of the equations in simulations, but in 
controls work the objective has usually been to ob- 
tain linearized EOMs for the application of linear sys- 
tem theory. Satisfactory results have been obtained 
for single-helicopter systems by this approach, but re- 
sults for dual lift have been limited to hover, where 
the second-order dynamics can be neglected. In ad- 
dition, modern control design methods are available 
for aeronautical systems with significant aerodynamic 
or dynamic nonlinearities based on partial inversion of 
the nonlinear EOMs. A segment of the dual-lift control 
literature seeks to apply these methods (refs. 17-21). 
For this work, a suitable nonlinear model is needed 
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which is both analytically tractable and computation- 
ally efficient. 

An alternative, ad-hoc analytical approach for in- 
elastic cables is based on manipulation of the Newton- 
Euler equations to eliminate cable tensions in a subset 
of these equations. This subset is combined with kine- 
matic constraint equations to obtain 6n equations for 
the linear and angular accelerations of the system’s 
rigid bodies, where n is the number of bodies. Re- 
sults have been obtained for single- helicopter systems 
(ref. 16) and, recently, for dual lift (twin lift), (refs. 19 
and 20). This approach requires the inverse of a larger 
system matrix than is required for generalized coordi- 
nates, and the choice of coordinates is restricted. 

Evidently, the dynamical equations for both elastic 
and inelastic suspensions are of interest. The existing 
formulations for these two suspension models differ sig- 
nificantly in form and in analytical and computational 
requirements, and difficulties of derivation or compu- 
tation associated with existing methods for inelastic 
cables become significant for dual-lift and multilift 
systems. 

In the present work, slung-load dynamics are treated 
systematically. For this purpose, many formulations 
of the equations for dynamical systems can be ap- 
plied to both elastic and inelastic suspensions. The 
Newton-Euler equations are suited to slung-load sys- 
tems and are used herein. Simulation equations are de- 
rived by applying the Newton-Euler equations to each 
rigid body, defining configuration vectors for the 77 - 
body system, and introducing generalized velocity co- 
ordinates. The result for inelastic cables is obtained by 
applying d’Alembert’s principle. Two formulations are 
obtained, one for elastic and one for inelastic suspen- 
sions. These formulations generalize the previous case- 
specific results to general slung-load systems and any 
set of generalized coordinates. The method is similar 
to that described in references 22-24 for general multi- 
body systems. The result for inelastic cables requires 
the inverse of a d x d matrix for which an analytical 
inverse is unknown, where d is the number of degrees 
of freedom (DOF). Slung-load systems with inelastic 
cables have only a few constraints, so d is relatively 
large, near 6?r in all cases, where n is the number of 
rigid bodies. Values of d are listed in figure 1; d = 20 
for the dual-lift system 10. 

A second pair of formulations is obtained by select- 
ing the generalized coordinates to represent the con- 
straints found in slung-load systems. The general- 
ized coordinates for the elastic system are selected to 


be partitioned into d coordinates that represent the 
configuration motion with invariant cable lengths and 
c coordinates that define the motion caused by cable 
stretching, where c is the number of constraints im- 
posed by the inelastic suspension. This can usually 
be done by including appropriate cable velocity co- 
ordinates in the generalized coordinates. The result 
for inelastic cables is obtained by nulling the cable- 
stretching coordinates to obtain an equation for the 
resultant constraint forces on the configuration’s rigid 
bodies. The solution for this equation requires the in- 
version of a relatively small c x c matrix (c = 4 for 
dual-lift system 10 of fig. 1), and the resultant force 
appears in the EOMs as an additive force, just as in 
the elastic-cable equations. 

The second formulation has reduced the computa- 
tional penalty relative to the elastic-cable equations, 
can be integrated with the elastic-cable equations in 
one simulation, and can be expanded nearly com- 
pletely in terms of natural vectors and matrices to 
obtain compact, tractable formulations. It is read- 
ily applied to the multilift systems. Results are given 
herein for the three dual-lift systems of figure 1 (sys- 
tems 8-10), and these can be integrated into a single 
simulation. Additional results are included for some 
single-helicopter systems to permit comparison with 
past work, and for the multilift system (system 11) of 
figure 1 extended to any number of helicopters. The 
derivation and results are sufficiently brief so that hand 
analysis and computer programming are practical. In 
addition, degenerate-body (point masses, rigid rods) 
approximations are given for general slung-load sys- 
tems and for two multilift cases. Last, linearized equa- 
tions of motion are formulated for general slung-load 
systems, and for dual-lift system 10 of figure 1. 

Equations of Motion for Multibody Systems 

A large body of literature on the dynamics of multi- 
body systems has accumulated since the early 1960s 
in response to the increasing importance of multibody 
dynamics in the design of, for example, spacecraft, ma- 
chines, robotic arms, and human motion models, and 
the relevance of the approaches used and of the results 
obtained to multilift helicopter systems is of interest. 
The principal aim in the literature has been to develop 
general-purpose computer programs to provide EOMs 
from a minimal amount of user input data defining 
the multibody system. This aim is motivated by the 
impracticality or excessive labor and unfavorable error 
probabilities of using hand derivation in most working 
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circumstances in these applications. Theory and anal- 
ysis are given, for example, in references 24-26, and 
surveys of computer programs that encode EOMs for 
a general system or that generate and compile sym- 
bolic case-specific code from user inputs are provided 
in references 27 and 28. 

Slung-load systems differ from the applications dis- 
cussed above. First, slung-load systems with inelastic 
cables have only a few constraints relative to the num- 
ber of DOFs (c << d, 6n) and can be equally well 
represented as unconstrained (elastic cables). The ap- 
plications cited above are^all highly constrained, with 
relatively few DOFs. For example, spacecraft and 
robotic arms are commonly represented as n rigid bod- 
ies with fixed orbit or base connected by n — 1 joints 
which permit one DOF of relative rotational motion; 
hence, d = d b + n — 1, where d 5 is the number of base 
body DOFs. Consequently, the conventional formula- 
tion in which a d x d matrix is formed and inverted is 
computationally more efficient for these applications 
than one containing a c x c matrix, but the converse 
is true for slung loads. Moreover, in most references 
the interbody connections are modeled as joints, and 
then convenient generalized coordinates are predefined 
according to the joint model. These have little appli- 
cability to slung-load suspensions. Gimbal-type joints 
can represent one to three DOFs of relative rotational 
motion, whereas inelastic suspensions allow three to 
five DOFs of relative motion, and elastic suspensions 
impose no constraints. 

One code-generating program, NEWEIJL, and its 
underlying formulation (given in refs. 22-24), does not 
specialize the interbody connection and can be applied 
to slung loads. The equations from this program are 
in the conventional form which requires the inverse of 
a d x d matrix, and the results are obtained with all 
terms expanded to their scalar elements. 

An alternative computer-based approach is to use 
symbolic computations to carry out routine analyti- 
cal steps — for example, the energy derivatives in La- 
grange’s equations. Some general possibilities of ap- 
plying MACSYMA (ref. 29) for this purpose are dis- 
cussed in reference 30 and this approach is used in ref- 
erence 7 to obtain linearized hover equations for dual- 
lift systems from Lagrange’s equations. 

Symbolic computation has many practical advan- 
tages, especially that of error-free derivations, and it 
has been applied where possible in the present work. 


Its ability to generate useful scalar equations for dual- 
lift systems from traditional methods has been lim- 
ited by the explosion in the number of scalar dynamic 
terms. In the present efforts, it appeared necessary 
to seek more efficient analytical methods. A new for- 
mulation and a new applications technique are pre- 
sented which improve the computational efficiency of 
the inelastic-suspension equations over previous forms, 
and which render hand derivation, analysis, and pro- 
gramming feasible for the previously difficult dual-lift 
and multilift systems. 

2. SYSTEM DESCRIPTION 

The systems of interest consist of one or more heli- 
copters that support a load (or more than one load, in 
some instances) by means of a suspension. For typical 
slung loads and nominal trajectories, the total load to 
be supported by the helicopters due to load weight, 
acceleration, and aerodynamics, is dominated by the 
load weight. The suspension consists of cables, usually 
of nylon webbing, and hooks, rings, isolator springs, 
spreader bars, and other hardware (ref. 31). Sus- 
pension designs with controllable geometry obtained 
by active cable winching and attachment- point move- 
ment have been proposed for load- attitude stabiliza- 
tion (refs. 32 and 33), and are included in the present 
formulation. 

These systems can be partitioned into n rigid bod- 
ies of non-negligible weight, (helicopters, 

load, spreader bar). The remainder of such sj^stems, 
referred to as the “suspension” hereinafter, consists 
of 771 straight-line links which support force only in 
the direction of the link (or only tension, in the case 
of cables) and have negligible mass and aerodynamic 
force compared to those of interest in trajectory con- 
trol. The number of suspension links m is listed in 
figure 1 for the systems shown there. In many cases, 
short links attached to the load are considered part 
of the load rigid body, because their elastic stretch is 
negligible and they are immobile relative to the load. 
The links can be modeled as inelastic, in which case 
c < m holonomic (position) constraints are imposed 
on the motion of the rigid bodies, and the system has 
d = 6n - c DOFs. Values of c and d are listed in 
figure 1. Alternatively, the links can be modeled as 
elastic, because of cable or isolator-spring stretching, 
in which case there are 6n DOFs. 
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An examination of cable and suspension elasticity 
and its effect on rigid-body motion is reported in ref- 
erence 34. Cable stretching under tension is usually 
modeled as that of an undamped spring with damp- 
ing supplied by the aerodynamic resistance of the at- 
tached bodies. Cables tend to be stiff, but the sus- 
pension design must avoid an upper bound which is 
set by resonance with the helicopter rotor frequency 
(about 4-5 Hz) where a divergent pilot-induced verti- 
cal bounce mode has been observed near hover. The 
net result is that natural frequencies of practical sus- 
pensions are about 2-2.5 Hz. This frequency is suf- 
ficiently high to be disjoint from the frequency range 
of interest in trajectory control (about 0.5 Hz). The 
corresponding mode is one of rapid and significant ca- 
ble tension variations, but with small stretching ex- 
cursions so that the rigid-body coordinates are nearly 
unaffected. 


3. EQUATIONS OF MOTION OF GENERAL 
SLUNG-LOAD SYSTEMS 

Configuration Vectors 

Physical vectors are referred to inertial or body-axis 
reference frames in the following discussion. Trans- 
lational motion and forces are given in inertial co- 
ordinates, and rotational motion and moments are 
given in body axes. The reference frame is indi- 
cated by a subscript , which is N for inertial space and 
i E {1, 2, . . . , n} for body axes of the zth body. Body- 
axis components of translational velocity and motion 
variables relative to a reference body are commonly 
used in slung-load simulations, and are readily intro- 
duced later when generalized coordinates for an appli- 
cation are selected. 


Both elastic- and inelastic-suspension models are 
of interest in trajectory control. For practical sus- 
pensions, simulations can employ the simpler, more 
general, and more computationally efficient nonlinear 
equations of the elastic model. If actual suspension 
stiffness were significantly greater, then difficulties of 
numerical stability and ill conditioning might arise in 
real-time digital simulation of the higher-frequency, 
low-amplitude cable-stretching motion. In control de- 
signs, practical suspensions can be approximated as in- 
elastic thereby eliminating feedback of states that have 
negligible influence on the rigid-body motion. If actual 
suspension stiffness were significantly lower, then the 
lower-frequency, higher- amplitude motions caused by 
cable stretching would be of interest in trajectory con- 
trol. In the present report, simulation equations for 
both elastic and inelastic suspensions are considered. 

The simulation of cable collapse is an application 
detail outside the scope of the present paper. If a ca- 
ble collapses, the resulting system is still a member 
of the class of interest and can be simulated. Practi- 
cal suspensions are designed and operated such that 
cable collapse does not occur except during large un- 
stable excursions from the nominal configuration, or 
during load pickup and dropoff. System 4 of figure 1 
(inverted-V suspension) is a special case in which one 
or two cables can collapse during small yawing motions 
if the cables are assumed inelastic, but this does not 
occur in cables with typical elastic properties. 


It is convenient to use configuration vectors that 
define the motion and forces of the n rigid bodies 
whose masses, inertias, c.g. translational motion, Eu- 
ler attitude-angle triplets, and angular velocities rel- 
ative to inertial space are denoted by (ml, 

VI*, al, wl),. . . ,(mn, Jn,Rn*, Vn*, an,cjn). The 
configuration vectors of position r and velocity v are 
defined as lists of the rigid body c.g. positions and 
Euler attitude angles, and the rigid-body c.g. transla- 
tional and angular velocities: 


/ ^l*A r \ 


Rn* N 

al 


\ an / 


(VI *N\ 


Vn* N 

ul i 


V um n ) 


(i) 


Let / be a list of the resultant forces and mo- 
ments applied to each rigid body, and let fg , /a, fc 
be corresponding lists of the gravitational forces, the 
aerodynamic and rotor forces, and the cable forces, 
respectively: 

f = f9 + fa + fc 


TT 
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( m\ gN \ 


/ FAIn \ 

fg = 

mn gN 
0 

/a = 

F An n 
MA h 


\ o ) 


\MAn n J 


/ FC\ n \ 


fc = 


FCtim 

MCh 


( 2 ) 


\MCn n / 

where FAi, MAi are the sum of the aerodynamic and 
rotor forces applied to the ith body and the sum of 
their moments about its c.g., respectively; and where 
FCi, MCi are similar force and moment sums due to 
all cables acting on Bi. The applied forces / depend 
on r, v, and the helicopter controls; any dependence of 
the aerodynamics on the acceleration v is assumed to 
be negligible. 


Last, let /* be a list of the inertia reactions of the 
n bodies: 



{ ml V1* N \ 


0 1 

r = - 

mn Vn* N 
J 1 cLTi 

- 

0 

S(u>h)Jl 


V Jn um n ) 


\ J TVLOTlyi / 


= -D v-X (3) 

The term f* is abbreviated to the form -D v-X, 
where D is block-diagonal with masses and inertia ma- 
trices along the diagonal, v is the configuration accel- 
eration, and X contains Coriolis terms due to the use 
of body-axis components of rotational motion. 


Kinematics 

The systems of interest consist of n rigid bodies con- 
nected by m cables which impose c constraints on the 
motion of the rigid bodies; they have d = 6 n - c 
DOFs, where c = 0 for elastic cables. These 

are holonomic systems; that is, the constraints im- 
posed by an inelastic suspension can be given as func- 
tions of position only. These constraints are usually 


time-invariant, but in the special case of active ca- 
ble winching or attachment-point movement they have 
explicit dependence on time. To accommodate this 
case, the kinematic model below includes parameters 
p = (pi,...,p 7r ) T , which can represent the control- 
lable geometric parameters of such a suspension or 
other known time-varying parameters convenient to 
the kinematic model in particular applications. 


For holonomic systems with d DOFs there exist 
d generalized position coordinates, q = {qu 9d) > 
which suffice to locate all points in the system and also 
the configuration position 

r = r(q,p ) (4) 

and d generalized velocity coordinates which suffice to 
define all inertial velocities of the system: 

u = U{q,p)q (5) 

The configuration velocity is related to u by a linear 
expression of the form 

v = A(q, p)u + B(q, p)p (6) 

Here, U is a nonsingular d x d matrix; it can be unity, 
but velocity coordinates different from q are commonly 
useful in applications. Note that v is asserted in equa- 
tion (6) to be linear in u.p. This follows from the 
usual linear relationship v(r) from rigid-body kinemat- 
ics (appendix A), 

«'» = *•» \ i = i „ 

t oii = WU{ai) ai J 

and equations (4) and (5). The term .4 is a 6n x d 
matrix. For inelastic suspensions, A expresses the con- 
straints by confining that part of the configuration ve- 
locity due to u to the d-dimensional linear vector sub- 
space defined by the columns of A. If p = 0, then 
this subspace is tangent to the configuration trajec- 
tory r(t). 

Equations of Motion 

The Newton-Euler equations for each rigid body are 
as follows: 

Fin — mi Vi* N = 0 
Mii — Ji c uii — S(udi) Ji — 0 

where Fi, Mi are the total applied forces and moments 
about the c.g. of body Bi. The same equations listed 
for all n bodies are 

fc+ fg + fa + f* =0 (7) 
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To obtain the simulation equations, differentiate equa- 
tion (6) with respect to time, introduce the result in f* 
(eq. (3)), premultiply equation (7) by A T , and solve 
for u : 

* = [A T DA]~ x A T {fo + fc 1 (8) 

where 

fo* fg + fa - DAu-X- D(B p-h Bp) 

The configuration vector fo denotes the combined ex- 
ternal forces, second-order velocity effects resulting 
fiom the choice of coordinates u, v, and the inertia 
reaction of the configuration to p(t). 

If the cables are elastic, then A is a nonsingular 
6n x 6n matrix and 

u = A~ 1 D~ 1 [fo + fc] (9a) 

In the case that p = 0 and we choose u = v, then 
^ I and the result is identical to the Newton- Euler 
equations applied to each body: 

V = D~ l [fg + fa- X + fc] (96) 

Equation (9b) for the rigid-body velocity coordinates 
can be applied to any configuration without further 
analysis except as needed to express fc. Equation (9a) 
generalizes this case to allow a choice of generalized 
velocity coordinates; for example, the use of cable- 
velocity coordinates in u provides a convenient and 
well-conditioned calculation of cable lengths and di- 
rections for use in the calculation of the cable forces 
fc. It should be noted here that the matrix A~~ x in 
equation (9a) represents the kinematic relation u(v) 
and can be given analytically from the kinematics as 
readily as the matrix A representing the reverse rela- 
tion v(u). Therefore, it is unnecessary to perform a 
matrix inversion to obtain A" 1 . 

The constraint force fc in equations (9) can be given 
as a sum of forces and moments applied by the suspen- 
sion at each attachment point: 

M 

fc=Y.hj TCj 
1 

where j enumerates every attachment of a cable to a 
rigid body, hj is a configuration vector defined in the 
next section, and TCj is the cable tension, which is 
given by the spring model of the cable as 

TCj = max{ 0, Kj{lj-toj)+cj tj} j = 1, 2, . . . , M 

where {(oj, Kj, cj) are the unloaded cable length and 
cable spring and damping constants, respectively Ca- 
ble damping, cj ^ 0, is introduced in reference 15, 


but otherwise has been neglected in simulations with 
elastic-cable models. For cables made of nylon web- 
bing, the spring rate K depends on loading and cable 
length (ref. 34). This dependency is readily included in 
the calculation of cable tensions from the above equa- 
tions. If accurate simulation of the small motions that 
result from suspension stretching is unimportant, then 
a simplified spring model of the cable elasticity that is 
consistent with typical natural modes will suffice. The 
spring constants K, c are related to the cable’s natu- 
ral mode parameters by a ~ Kg/F , 2£u/ n = cg/F , 
where F is the load supported by the cable. In simula- 
tions with elastic cables, initial values of cable tension 
and stretch can be calculated in cases with nonredun- 
dant suspensions from the solution for the constraint 
forces fc for inelastic cables given below. Alterna- 
tively, the configuration can be allowed to settle from 
an approximate initial arrangement, possibly with the 
aid of cable damping as a settling device (ref. 15). 

Simulation Equations for Inelastic Suspensions 
Using d’Alembert’s Principle 

If the cables are inelastic then the cable forces 
fc drop out of equation (8) {A T fc = 0) in ac- 
coi dance with d’Alembert’s principle for constrained 
holonomic systems; the constraint forces do no virtual 
work. This result is shown as follows. First, enumer- 
ate the cable attachments 1,2, .. . ,M at attachment 
points Rl,R2, . . . ,RM on their corresponding bod- 
ies Bi , i = i( 1), i( 2), . . . , i(M). One or more cables are 
attached at an attachment point and every such at- 
tachment of a cable to a body is numbered. The con- 
straint force on the configuration is then (see fig. 2(a)) 

M 

fc = J2 TC J ( 10 ) 

i 

where TCj is the cable tension. The jth cable at- 
tachment at point. Rj on body i(j) applies a force and 
moment to Bi{j) given by 

FCij = kcj TCj 
MCij = (Ri*j x kcj) TCj 

where kcj is the cable direction outward from the body 
and Ri*j is the moment arm of the attachment point 
about the c.g. (Rj — Ri*). Thus, hj in equation (10) is 
a configuration vector whose nonzero elements are kcj 
and (Ri*jxkcj), corresponding to the constraint force 
and moment on Bi(j) due to the jth cable attachment, 
per unit tension. 

Second, enumerate the cables and links which con- 
stitute the suspension, Cl, C2, . . . , Cm. Each end of 
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(b) Cable connecting 
two rigid bodies 

Figure 2. Suspension forces. 


a cable or link is attached to either a rigid body or to 
another link. In the restricted case, in which all cables 
are attached at both ends to rigid bodies, all terms in 
equation (10) can be combined in pairs, with each pair 
corresponding to a single cable: 

m 

/c = ^FjTCj (11) 

3 = 1 

where the nonzero elements of Hj are 
kcj, — kcj, (Rj — Ri*) x kcj, — (R j 7 — RF*) x kcj, 
corresponding to the forces and moments on bodies 
Bi , Bi' to which the cable Cj is attached at points 
Rj, Rj' (see fig. 2(b)). Here, kcj is the direction of 
the line segment (Rj' — Rj)- From the definitions of 
Hj, v it follows that v T Hj combines the rigid-body 
c.g. and angular velocities to give the difference in 
cable endpoint velocities along the cable direction: 

v T Hj = Vi* • kcj — Vi * • kcj + u>i • (Ri* j x kcj) 
- u)i l • {Ri*f x kcj ) 

= (Vj - Vj') • kcj = -tj 

That is, it is the cable stretch rate. Recalling equa- 
tion (6), u T A T Hj is the cable stretch rate due to u. 
Since all cable lengths are invariant for arbitrary u in 
an inelastic system, then 

A T Hj = 0 j = 1, 2, . . . , m 


and 

A T fc = 0 

For inelastic suspensions, all systems in figure 1 except 
system 6 can be partitioned into rigid bodies such that 
every cable connects two rigid bodies. 

For more general suspensions in which links are con- 
nected at both ends to a rigid body or to another link 
(fig. 2(c)), fc is given by equation (10), and v T fc is a 
linear combination of the attachment-point velocities. 
To this combination, first add and subtract the veloc- 
ities of the cable interconnection points in the cable 
directions, and then apply the force-balance condition 
to the linear combinations of cable forces at these in- 
terconnections. The result is 
M 

v T fc = ^ Vj . kcj TCj 

m 

= H(Vj-Vj> kcj TCj 

3 = 1 

771 

= -]TfjTCj (12) 

3 - 1 

where the second and third sums are taken over all 
cables and {tj} are cable lengths. Consequently, if 
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the suspension is inelastic then all cable lengths are 
invariant for arbitrary u, and 

A T fc = 0 (13) 

Equation (13) indicates that the cable force fc is or- 
thogonal to all the columns of A and therefore to the 
configuration velocity v y when p = 0. This result also 
expresses d'Alembert’s principle for constrained holo- 
nomic systems (the constraint forces fc do no virtual 
work) and Jourdain’s principle for nonholonomic sys- 
tems (the constraint forces have no virtual power): 

fc T Sr = 0 

and 

fc T Sv = 0 

where the virtual position and velocity in this context 
are 

Sr = [VT v] 6q = AU Sq 

and 

Sv = [Vj u] Su = A 6u 

and the constraint forces to which these principles ap- 
ply in slung-load systems are shown above to be the re- 
sultant cable forces if the suspension is inelastic. Equa- 
tion (8) for inelastic suspensions is now 

u = [A T DA]~ l A T fo (14) 

Here, [A? DA\ is a positive definite d x d matrix; D 
contains the system mass-inertia parameters; A con- 
tains the system geometry and constraints; and fo (de- 
fined by eq. (8)) contains the term DA u -f X , which 
is second-order in velocity coordinates, and the iner- 
tia reactions to p. It is unnecessary to calculate cable 
tensions in this result. 

The principal difficulty in implementing equa- 
tion (14) is the required inversion of [A T DA]. A gen- 
eral analytical inverse is unknown. For slung loads, 
d is relatively large (near 6n) and the inversion is 
therefore computationally more costly than it is for 
the elastic-cable case (eq. (9)), which requires no in- 
version. Some analysis of its numerical conditioning 
over the flight envelope is needed before inversion can 
be confidently implemented. In addition, the matrix 
inversion obstructs further expansion and rationaliza- 
tion in terms of three-dimensional dynamics whereas 
equation (9) can be fully expanded in applications in 
terms of natural vectors and matrices to obtain an- 
alytically the sensitivities of the elements of ii to the 
forces and moments. A Gauss-Jordan elimination pro- 
cedure to reduce the matrix to a block-diagonal form 


was considered in initial efforts to reduce the inver- 
sion to the inverses of smaller blocks given in terms 
of three-dimensional vectors. The procedure was ap- 
plied in each case after defining the coordinates u and 
identifying some subdivision of u into natural vectors 
or groups of coordinates with which to associate the 
blocks of the diagonalized form. For typical choices 
of u, these blocks usually represented the rigid-body 
inertia matrices modified by the constraint moments 
(effective inertia). 

Analytical results were obtained for systems 1, 2, 
7, and 8 of figure 1, but the labor required increased 
rapidly with the number of DOFs, and the procedure 
appeared impractical to do by hand for dual-lift sys- 
tem 10 of figure 1. No generalization of the effective 
inertia interpretation of the blocks or of the simplifica- 
tion required at each step of the reduction was found, 
so the procedure was ad hoc and specific to each sys- 
tem and choice of coordinates, and could not be given 
effectively as a MACSYMA algorithm for analytical 
reduction. 

A second form of the simulation equations, which 
requires the inverse of a much smaller cxc matrix, was 
obtained after restricting the generalized coordinates 
appropriately for slung-load systems; it is described 
next. The result can be applied with much less labor 
and likelihood of error than the procedure considered 
above. 

Equation (14), or equivalent forms obtained by 
means of Lagrange’s equations, is implemented in 
many of the general multibody programs, where d is 
relatively small in typical applications in the multi- 
body literature. Numerical inversion of positive defi- 
nite matrices can be carried out by various algorithms 
based on the Cholesky (square-root) decomposition 
(ref. 35) previously applied in the field of estimation. 
These are computationally efficient and resist numer- 
ical instabilities. The conditioning of the coefficient 
matrix \A T DA) depends on A in equation (14) or, 
equivalently, on the choice of coordinates u. In the 
multibody programs these coordinates are often pres- 
elected based on the joint model of interbody connec- 
tions, and these appropriate coordinates tend to result 
in well-conditioned coefficient matrices in typical ap- 
plications. In the present work, the interbody connec- 
tions are suspensions composed of cables; it is left to 
the applications phase to determine in each case what 
constraints are imposed and what choice of coordinates 
best represents the constrained system motion. 
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Simulation Equations for Inelastic Suspension: 
Explicit Constraint Method 

The approach in this section is to select the gen- 
eralized velocity coordinates for the elastic system in 
equation (9) so as to separate the motion due to ca- 
ble stretching from the remaining motion by means of 
invariant cable lengths, and then to study the results 
when the stretching motion is nulled. A solution for 
the constraint force and new simulation equations for 
the inelastic system are obtained. 


where fo is as defined with equation (8) except that 
Al ul replaces A u. As in equation (14), fc drops 
out (Al T fc — 0 in view of eq. (13)) and the result 
differs from the inelastic system equations (eq. (14)) 
only in the presence of the configuration acceleration 
due to elastic stretching. The effect of the stretching 
motion on ul depends on the cable spring constants: 
for a fixed disturbance, the extremes of A decline with 
increasing spring stiffness, whereas the extremes of the 
term in A remain fixed in magnitude, although A (t) 
increases in frequency (ref. 15). 


First, let the generalized velocity coordinates of the 
elastic system u be composed of 6n — c coordinates, ul, 
of the system with invariant cable lengths, and c length 
rates, A, which suffice to define the motion resulting 
from cable stretching. In general, if the c independent 
position constraints imposed by an inelastic suspen- 
sion are given by {Al(r) = 0, . . . , Ac(r) = 0} then A 
can be taken as (Al, . . . , A c) T . For slung-load systems, 
A can usually be taken as the cable lengths and the 
complete set of coordinates (ul, A) can be taken as the 
c.g. velocity of a reference body, the angular velocities 
of all rigid bodies, and the cable angular velocities and 
stretching rates. 

Next, substitute the partitioned u in equations (6) 
and (9): 


The constraint force can be given in terms of c in- 
dependent parameters; that is, it can be given in the 
form 

( Sl \ 




s2 


c J 


= H s 


where s is arbitrary, {Hj} are configuration vectors, 
and rank [H] — c. From equation (13), Al T H ~ 0, 
and from the construction of A -1 , Al T A = 0. There- 
fore, the columns of H and A are both bases of the 
same linear vector space, and A can be used to define 
the constraint force: 


fc = A s 


(17) 


ul 


— Al ul 4" L A B p 

(15) 

'AI1 t \ 


D~\fo + fc) 

(16) 

, A T / 



where fo is defined in equation (8); Al, L are the 
6n — c and c columns of A which, respectively, define 
the contributions of ul, A to u; and Al 1 T , A T are the 
6 n — c and c rows of A” 1 which, respectively, define the 
relations ul(u) and A(u). As noted earlier, A -1 can be 
obtained without matrix inversion since it defines the 
relation u(v) for the elastic system and can be given 
from the kinematics as readily as A, which defines the 
relation u(u). 


Equation (16) gives the simulation equation for sys- 
tems with elastic suspensions in terms of the coordi- 
nates (ul, A), where ul leaves the cable lengths in- 
variant. The influence of cable-stretching motion on 
ul can be viewed by entering the derivative of the 
partitioned generalized velocity coordinates given by 
equation (15) in /*, and rederiving equation (8). The 
first 6n — c equations can be arranged as 

ul = [A1 T D Al]" 1 Al T [fo - D(LX + LA)] 


where A r is the Jacobian [Vj A], and s has units of 
force, if the coordinates A are lengths. 

For inelastic suspensions, A = 0, and equation (16) 
gives c equations from which to calculate the con- 
straint forces: 


0 = A t D 1 (fo + fc) (18) 

Introduce equation (17) into equation (18) and solve 
for s: 

s = -[A T D" 1 A]- 1 A T D- 1 /o (19) 

Further, A can be replaced in equations (17) and (19) 
by any other convenient basis, of the space con- 
taining fc . For example, in the special case that every 
cable connects two bodies and c = m, then s can be 
taken as the c cable tensions with the basis vectors 
H 1, . . . , LTc, as defined by equation (11) above. 

The cable tensions are related to s, and can be 
uniquely determined from s if the suspension is not 
redundant. However, the constraint force applied to 
the configuration, /c, can always be calculated. A 
suspension can be separated into disjoint sets of inter- 
connected links. Each such set imposes one holonomic 
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constraint; if the number of sets is c, then all cable 
tensions can be uniquely determined from s, but if it 
exceeds c, then the constraints are redundantly im- 
posed. In the special case that all cables connect two 
rigid bodies, then each cable is a disjoint set and all 
cable tensions can be found only if c = m. 

Finally, the simulation equations are 

u — A~ 1 D~ l (fo + A s) (20) 

or, for the inelastic DOFs; 

ill = AIl T D- l [I - h\i\ T D- l i\]~ l A T D- l ]fo (21) 

This result has several advantages over equation (14) 
for the slung-load systems. First, the leading coeffi- 
cient matrix, A~ l or All , is known analytically. Sec- 
ond, equations (19) and (20) require the inversion of a 
c x c matrix, [A T Z> _1 A], where c — 1 for system 1 and 
c = 4 for dual-lift system 10 (fig. 1), compared with 
the required inversion of 11 x 11 and 20 x 20 matrices, 
respectively, for these two systems when equation (14) 
is used. In many cases, the cable tensions can also be 
generated by equation (19). 

Computational efficiency in calculating the dynam- 
ics can be compared among the formulations of the 
simulation equations defined above by counting the 
number of multiplications required to generate u, given 
A , A~ 1 , D, D ~ 1 , / 0 . These are shown in table 1 for sys- 
tems 1 and 10 of figure 1. 

In table 1, g{) refers to the number of multiplica- 
tions and divisions required for the Cholesky inver- 
sion, which increases with the square of the matrix size 
(ref. 35). The general formulas (derivations omitted) 
include savings gained from generic zeros, ones, and 
matrix symmetries. The number of coordinates u,v 
that are identical is represented by i in these formulas, 
and accounts for the generic zeros and ones. In most 
cases, u can be selected such that i = 3n+3. Additions 
are omitted from the operations count for brevity, but 
this omission does not affect the conclusions. Equa- 
tion (14) is representative of previous single-case for- 
mulations of the slung load dynamics with inelastic 


cables, as well as those in the multibody literature 
(e.g., refs. 24 and 36). As shown, equations (19) and 
(20) provide a significant reduction in the computa- 
tional requirements to represent these dynamics, and 
a much reduced penalty compared to the elastic cable 
formulation (eq. (9)). The same conclusion applies to 
all systems of figure 1 . 

Computational requirements are of interest in real- 
time simulation and control. In simulations, the time 
required to compute the dynamical terms are of in- 
terest if these are a significant fraction of the cycle 
time. In control systems, linear control laws use coef- 
ficient matrices that can be calculated offline, and are 
nearly unaffected by the underlying formulation of the 
nonlinear EOMs. However, modern control-law design 
methods for aeronautical systems (inverse model, out- 
put linearization) treat systems with significant aero- 
dynamic or dynamic nonlinearities by partially invert- 
ing the simulation equations. It is beyond the scope 
of this report to consider such control laws in detail, 
but it is expected that the formulation of the dynam- 
ics in the EOMs will affect the computational require- 
ments for such control laws significantly for slung-load 
systems. 

Simulation Equations for Inelastic Suspensions 
from Lagrange’s Equations 

Lagrange’s equations have been the principal ana- 
lytical approach in slung-load controls studies, where 
the object is usually to obtain linearized EOMs. The 
EOMs for general slung-load systems with inelastic 
suspension are derived here from Lagrange’s equations, 
and it will be shown that they have the same form as 
equation (14). 

Lagrange’s equations for general holonomic systems 
with d DOFs are 

j t Vg KE - V, KE = Q (22) 

where q = (<?i, . - - ,qd) T are generalized position co- 
ordinates; Q = (Qi, . . . , Qd) T are the generalized 
forces, including conservative forces; and KE is the 


Table 1. Multiplications required to calculate ii 


Equation 

Multiplies 

System 1 

System 10 

General formula 

(9) 

(14) 

(19) and (20) 

60 

616 + <?(11 2 ) 
119 

264 

3430 + 5 (20 2 ) 
868 + g( 4 2 ) 

6n(6n — i + 2) 

d((6n i) ( -^ + 12n + d) + g(d 2 ) 

d(6n - i) + 3 n(c 2 + 9c + 4) + c 2 + g{c 2 ) 
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kinetic energy. For slung-load systems, KE , Q can be 
given as 

KE = (1/2) v t Dv = (1/2) q T A T DAq 

Q = A T (fg + fa ) (23) 

where A defines v(q) as 


v — Aq 


and v,fg,fa are defined in equations (1) and (2). 
The controllable parameters, p, in equation (6) have 
been omitted from the velocity relation and kinetic en- 
ergy for brevity, but they can be routinely included 
in KE and its derivatives. The generalized forces 
are obtained as coefficients of the generalized coor- 
dinate displacements in the virtual work, 8 W = f T 
6r = f T A8 q. Here, / is the resultant force on the 
configuration, / = fg + fa + fc. For inelastic cables, 
fc is orthogonal to the columns of A in accordance 
with d’Alembert’s principle, and drops out of Q. 


After carrying out the derivatives of the kinetic en- 
ergy in equation (22), the general slung-load equations 
are obtained as a second-order ordinary differential 
equation: 

M(q)q + k(q,q) = Q 
M = A t D A 

k = A T D Aq + (A- G) t D Aq (24) 


where 

G = [V^ (A(q) g)] 




dA . 

A Q* 
oqi 


All terms in k are second order in the velocity coordi- 
nates of q. Equation (14) can be placed in the same 
form as equation (24) for comparison: 


M(q)q + k(q , q) = A T {fg + fa) 


M — A T D A 

k = A T {DAq + X) = A T DA q + A T X (25) 

where u has been taken to be identical to q for the 
comparison. 


The equations of motion, (24) and (25), from La- 
grange’s equations and from equation (14), respec- 
tively, have the same form. The inverse of a d x d 
matrix, M , to obtain q is required in both cases. 


position coordinates q, forming KE(q,q ), and gener- 
ating the latter’s derivatives routinely (refs. 10, 15, 
and 7). Contact with the rigid-body velocities is lost 
in KE , and the terms M, A*, Q in equation (24) are 
obtained as structureless objects in the d-dimensional 
space of constrained motion. This has sufficed for 
studies of single helicopters, but appears to be imprac- 
tical for multilift systems. 

Equation (23) introduces the multibody structure of 
slung-load systems in the Lagrangian formulation, fol- 
lowing similar steps found in the multibody literature, 
to obtain equation (24) in terms of objects A, D, Aq , 
and / in the 6n-dimensional space of configuration mo- 
tion. Equations (24) and (25) are identical except for 
differences in the analytical statement of the second- 
order velocity terms, k. In particular, the term A T X 
in equation (25) is given easily from the kinematics 
and from equation (3), but its numerical equivalent 
in equation (24), (A - G) T D Aq , requires numerous 
derivatives of A. However, k can be neglected in the 
derivation of linearized equations for windless hover 
because its gradients with respect to q,q are zero for 
practical choices of coordinates, q. 

Additional differences between equations (14) and 
(24) occur in all terms when the velocity coordinates 
u are different from q. In the applications work, these 
terms can be generated more easily from equation (14) 
than from equation (24) owing to the simpler kine- 
matics A, usually obtained using generalized velocities, 
and the simpler expression for the nonlinear term, A;. 

Simulation Equations for Inelastic Suspensions 
Using Rigid Body Velocity Coordinates 

Simulation equations for the rigid-body velocity co- 
ordinates (eq. (1)) have been derived assuming in- 
elastic suspension for several slung-load systems in ref- 
erence 16 and for dual lift in reference 20. These spe- 
cialized derivations are obtained by extracting a sub- 
set of the Newton-Euler equations with the suspension 
forces eliminated and appending derivatives of the con- 
straint equations to obtain 6n equations for the rigid- 
body accelerations v. This procedure can be formally 
extended to general slung-load systems with inelastic 
suspension by using the partitioned coordinates previ- 
ously defined in equation (15). 

As before, partition the generalized velocity coordi- 
nates of the elastic system u into 6n — c coordinates 
ul of the system with invariant cable lengths and c 
length rates A with suffice to define the motion due 


In the slung-load literature, Lagrange’s equa- 
tions (22) have been applied by defining generalized 
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to cable stretching. Then the kinematics can again be 
partitioned as 


v = An — Al ul + LX 



Controllable parameters p have been omitted from the 
kinematics for brevity but can be included routinely 
in the derivation. Here, the columns of A and the 
rows of A” 1 are partitioned into [Al,L] and [All, A], 
The term Al is a basis of the space perpendicular to 
the suspension forces fc and A is a basis of the space 
containing these forces. The suspension forces can be 
expressed as in equation (17) in terms of A and c in- 
dependent constraint force parameters s . Introduce 
this into the Newton- Euler equation (9b) and premul- 
tiply the result by A T to obtain the dynamic equations 
(note that Al T /c = 0 and L T A = I): 

Al T Dv = Al T {fg + fa - X) (27) 

s = L t (Dv -fg-fa + X) (28) 

A derivative of equation (26) provides c addi- 
tional kinematic constraint equations for inelastic 
suspensions: 

A = A T y + A t v = o (29) 

Together, equations (27) and (29) comprise 6N scalar 
equations for the rigid-body accelerations in which the 
suspension forces have been eliminated. These yield 
the simulation equations in the form 

v ~ M(r)~ l z 

where 

/ A1 T D \ 

M -( r J <“> 

/ Al T (fg + fa- X) \ 

\ -A T v ) 

Equation (30) requires the inversion of the 6n x 6 n 
matrix, M , composed of 6n — c columns from D A and 
c rows from A -1 . General results for its rank or for 
efficient inversion algorithms are not available. The 
second-order velocity terms occur in X,h T v. The 
use of the rigid-body velocity coordinates results in 
calculation of relative motions from small differences 
of large numbers in a simulation. 


The suspension force parameters, s, can be com- 
puted from equation (28) after computing accelera- 
tions. Alternatively, equation (19) gives s in terms 
of the applied forces and the velocity coordinates, but 
requires more computation than equation (28). In gen- 
eral s suffices to define the resultant suspension forces 
applied to each rigid body of the configuration. Cable 
tensions are related to s by fc~ As = Ht , where t 
is a list of cable tensions and H contains appropriate 
configuration vectors defining the force and moment 
action of each cable on each rigid body. If the number 
of cables equals the number of constraints, then the 
cable tensions are readily obtained from this relation. 

4. APPLICATIONS 


The object of this section is to demonstrate and ap- 
ply the methods of this report to slung-load systems 
of practical interest, particularly the dual-lift systems 
for which general simulation equations were not pre- 
viously available. Results are given for systems 1, 2, 
and 7 and for all three dual-lift systems shown in fig- 
ure 1, as well as for multilift system 11 extended to 
any number of helicopters. These results are given as 
summaries of simulation equations in programmable 
form in the appendixes. 

General Procedure 

A general procedure for applying and implementing 
the present results for slung-load simulation equations 
is outlined in figure 3. Both elastic suspensions (eq. (9) 
or (16)) and inelastic suspensions (eqs. (19) and (20)) 
are included. The first task (fig. 3(a)) is to perform 
an analysis to (1) determine the constraints of the in- 
elastic suspension; (2) define 6n generalized velocity 
coordinates (ul,A) such that ul are d coordinates for 
the inelastic suspension and such that A are c coor- 
dinates which define system motion caused by cable 
stretching; and (3) obtain expressions for A, A” 1 , A 
from the kinematics and for B, B from the suspen- 
sion geometry, if it is controlled. These items are case 
specific. 

The selection of appropriate generalized coordinates 
u is case specific, but several features were used repeat- 
edly in the applications to maintain simplicity of the 
kinematics and equations. The coordinates u can be 
selected to consist largely or entirely of natural vectors. 
If u contains rigid-body velocities identical to those in 
v , then the corresponding rows of A, A -1 , A are from 


n 
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c number of holonomic constraints on configuration for inelastic suspension 

u l 6 n-c generalized velocity coordinates for system with inelastic suspension 

A c generalized velocity coordinates defining suspension stretching motion 

u 6 n generalized velocity coordinates for system with elastic suspension 
u = (ul T ,\ T ) T 

A the Jacobian [V£ v] obtained from the kinematics 

A~ l the Jacobian [Vj u] obtained from the kinematics 

A the Jacobian [V„ A T ] from A -1 or any basis of the interaction force space 

B the Jacobian [VT ul (required only if P ± 0) 

(a) Quantities required a priori. 


fa = fa(r,v,8) 

f 0 = fg + fa- D Au- X - D(Bp 4- Bp) 

(6 = helicopter controls) 

hj muo:{0, Kj (£j - £oj)} 

f c — \ j=i 

[ -\[A T D~ l A}- l A T D- l fo 

(elastic cables) 
(inelastic cables) 

sf = D~ l (fo + fc) 


u — A~ l sf 


u = f u dt 


v = A u + B p or z = Z u u + Z p p 


r = r(r,p,v,p) 


II 

o 

II 

_ 


(b) Simulation equations. 


Figure 3. Procedure for simulation of slung-load systems. 
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(c) Simulation flow diagram. 
Figure 3. Concluded. 


the unit matrix or zero. If u contains cable veloci- 
ties, the corresponding rows of A,A~ l ,A contain only 
coordinate transformations and skew-symmetric ma- 
trices representing Coriolis velocities and Coriolis and 
centrifugal accelerations. In most of the applications 
discussed in this paper, u consists of the c.g. velocity 
of a reference body, the angular velocities of all bod- 
ies, and various cable velocities or their equivalents. A 
consideration in choosing u is that the relation v(u) 
should be nonsingular and well-conditioned over the 
domain of motions (v,r) of interest; that is, all con- 
figuration motion should be readily detected from a 
knowledge of it. Among the examples given, this was 
a factor only for system 7 (fig. 1). 


The second task is to implement the simulation 
equations (16) or (19) and (20) as given in the equa- 
tion summary (fig. 3(b)) and the information flow di- 
agram for the general ?t-body system (fig. 3(c)). The 
applied forces and moments due to weight, aerody- 
namics, and rotor arc combined in fo with the inertia 
coupling associated with the choice of coordinates u,v 
and the controllable geometry, if any. The configura- 
tion vectors fg , fa, X , and D were previously defined 
in terms of physical vectors and matrices of the three- 
dimensional rigid-body mechanics of the system’s con- 
stituent bodies, in equations (2) and (3). The aerody- 
namic term fa need not be considered in detail here, 
except to note that it can be assumed to be a function 
of the configuration kinematics and the acroctynamic 
controls, 6. These dependencies arc considered in more 
detail later in deriving the linearized EOMs. 

The forces and moments applied by the suspen- 
sion to the configuration fc arc calculated from equa- 
tions (10) or (11) (elastic cables) or equations (18) 


and (19) (inelastic cables). Both elastic and inelastic 
cable models can be accommodated by switching for- 
mulas for the interaction force. Expressions are given 
for both cases in all examples. For elastic cables, fc is 
calculated from cable tensions and from configuration 
vectors defining the force and moment action of each 
cable on the n rigid bodies. For inelastic cables, fc is 
calculated from any basis of the constraint force space 
A, and the corresponding independent constraint force 
parameters s. The Jacobian, [V v X T ], from A~ l is al- 
ways a basis. When every cable connects two bodies 
and when m = c, as in the dual-lift systems, then 
the constraint force parameters s can be taken as the 
cable tensions, and the coefficient matrix defined in 
equation (11) can be used for both models: 

fc = Hs 

(max{0, Kj(£j - toj )}, j = 1, . . . , c) T 
(elastic cables) 

-\H T D- 1 H J- 1 H T D~ } fo 

(inelastic cables) 

Subsequently, the configuration-specific forces and mo- 
ments sf, and, finally, ii are calculated. The remain- 
ing steps generate the kinematics u , v , r, r routinely or, 
more generally, these kinematic variables can be ex- 
panded to any set of velocity and position coordinates 
i, 2 of interest after defining z(z,p , u,p) from the kine- 
matics and geometry. For inelastic cables, the coordi- 
nates A and their equations can be eliminated since A is 
theoretically zero. Alternatively, all Gn equations can 
be retained to permit simulation of both cable models; 
for inelastic cables, the computed A indicates numeri- 
cal accuracy. 
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The general simulation equations shown in figure 3 
are given in terms of vectors in the 6n-dimensional 
configuration space. These are expressed in all results 
in terms of the natural vectors and matrices of the un- 
derlying three-dimensional rigid-body mechanics and 
geometry. General formulas used to obtain this form 
are given in appendix A. The required equations are 
relations from the classic vector theory of rigid-body 
dynamics in a programmable form. First, the skew- 
symmetric matrix is introduced to represent the cross 
products which occur in the Coriolis velocities, Cori- 
olis and centrifugal accelerations, and cable moments 
that pervade A, A -1 , Au, and /c. Second, reference 
frames are defined for all rigid bodies, and standard ex- 
pressions are given for the coordinate transformations, 
angular velocities, and transformation rates that occur 
in the results. All transformations are given in terms 
of the customary Euler-angle sequence of aeronautics. 
Third, the treatment of cable direction angles, cable 
axes, and cable velocities is outlined; Euler angles are 
again used, but only two angles are needed. 

Examples 

The literature describes a variety of existing and pro- 
posed slung-load systems (fig. 1). These are grouped, 
for this discussion, into (1) single helicopters with a 
single attachment point on the helicopter, (2) single 
helicopters with multiple attachment points, (3) dual- 
lift systems, and (4) multilift systems. Several suspen- 
sions are of interest in each category, including con- 
trollable geometry in systems 5 and 6 of figure 1. Re- 
sults are given for systems 1,2, and 7 and for all three 
dual-lift systems, as well as for a simple multilift sys- 
tem. Systems with controllable suspension geometry 
are not included in the present examples. 

The single-helicopter systems 1, 2, and 3 of figure 1, 
with a single attachment point on the helicopter, are 
considered in appendix B. Results are obtained for sys- 
tems 1 and 2 by using the above procedure. The 
suspensions for systems 1 and 2 impose three con- 
straints and one constraint, respectively, on the load 
relative motion when inelastic. Appropriate partition- 
ing coordinates were readily obtained by including the 
load relative velocity (system 1) or cable-axis compo- 
nents of the cable velocity (system 2) in the gener- 
alized coordinates. Differences from the earlier for- 
mulations of simulation equations include the use of 
relative-acceleration coordinates for the elastic suspen- 
sion case, the nonlinear EOMs for the inelastic suspen- 
sion with explicit calculation of cable tensions, and 


the integration of both elastic and inelastic suspension 
models into one equation set. 

Partitioning coordinates for single-cable-and-sling 
system 3 of figure 1 were not found. If the suspension 
is inelastic, this system is identical to system 2 and can 
be represented by the same equations. For the elastic 
suspension, the principal elastic elements are the sling 
legs, in which case no cable connects two rigid bodies. 
Here, the suspension geometry and cable tensions are 
determined from the force balance at the sling vertex. 
Coordinates are given in appendix B that simplify this 
calculation relative to the method found in the earlier 
literature. Kinematic relations needed to obtain the 
EOMs from equation (9a) are given, but these depend 
on suspension elasticity parameters, as well as on ordi- 
nary kinematics. A complete set of EOMs is omitted 
because the coordinates given are not of the type that 
is of interest in the present applications work. In this 
case the EOMs for an inelastic suspension are sim- 
pler than those for an elastic suspension because of 
the complexity in calculating interaction forces when 
each elastic link is not connected to two rigid bodies. 

The single-helicopter systems, 4-7 of figure 1, with 
multiple attachment points on the helicopter, are con- 
sidered in appendix C. These suspensions were de- 
veloped in the 1965-1975 period to stabilize difficult 
loads such as the standard cargo container (MIL- 
VAN). They provide yaw restraint to stabilize elon- 
gated loads in a minimum-drag heading, and pitch at- 
titude is restricted. The addition of active cable-length 
and attachment-point control to increase load-motion 
damping was also considered. 

Results are given in appendix C for the bifilar sus- 
pension (system 7, fig. 1), which imposes two con- 
straints on load relative-motion when inelastic. This 
system also approximates the inverted-Y suspension 
(system 6, fig. 1) when the bar is sufficiently close to 
the load. For inelastic suspension, several sets of load- 
suspension coordinates have been used in the earlier 
literature. The load-suspension geometry is examined 
in detail in appendix C, and several choices of coor- 
dinates, some of which are ill-conditioned and some 
well-conditioned, are identified. One set is selected 
and expanded to a set of partitioning coordinates, and 
a complete set of nonlinear EOMs for both elastic and 
inelastic-suspension is given. Differences from the ear- 
lier literature include the use of relative accelerations 
for the elastic suspension case, and the formulation 
of efficient equations for the inelastic suspension with 
cable tensions explicitly calculated. 
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The remaining suspensions in this group are analyti- 
cally difficult to represent as inelastic. In the inverted- 
Y suspension (with or without the spreader bar) no 
cable connects two rigid bodies, and the force balance 
at the suspension interconnections must be considered 
as in system 3, discussed above. Further, if the bar or 
interconnections are close to the helicopter, the sus- 
pension can be represented as an inverted V, and if 
close to the load, as bifilar. For inelastic cables, there 
can be one, two, or three constraints, depending on 
these details of suspension geometry. Finally, it ap- 
pears that representation of the inverted-V suspension 
as inelastic is both unrealistic and intractable. Simula- 
tion equations for the elastic inverted-V suspension are 
given in the literature, along with approximate equa- 
tions in which cable stretching motion is neglected. 
Although the existing analytical difficulties in simu- 
lating the multipoint suspensions of practical interest 
are not relieved by the methods presented in this re- 
port, recall that the principal object here is to treat 
the difficult multilift systems discussed next. 

Simulation equations for the three- and four-body 
dual-lift systems (8-10, fig. 1) are given in ap- 
pendix D. Appropriate partitioning coordinates for ap- 
plying equation (16) were readily found by using cable- 
axis components of the cable velocities or their equiv- 
alents. The results integrate elastic and inelastic sus- 
pension models; since c = rn for these systems, the 
suspension forces fc= Hs can be calculated by using 
the identical matrix of configuration vectors H for the 
two models. Further, all three systems can be simu- 
lated by a single equation set. System 8 is a three- 
body subsystem of system 10 obtained by deleting the 
load and bridle cables, and can be represented by a 
subset of the coordinates and equations for system 10, 
and system 9 is a simple specialization of system 8 
with coincident attachment points on the load. Gen- 
eral nonlinear equations for these systems with inelas- 
tic suspensions were not previously available owing to 
the complexity of representing the dynamics by the 
traditional methods found in the slung-load literature. 
The methods of this report render these dynamics ac- 
cessible and tractable, and provide a new formulation 
of the EOMs for inelastic suspensions that is efficient 
and that provides explicit calculation of the suspension 
forces. 

Simulation equations for the multilift system 11 of 
figure 1 extended to any number of helicopters are 
given in appendix E. Each helicopter is connected to 
the load by a single cable, which, if inelastic, imposes 
one constraint on the motion of a helicopter relative 


to the load. Equations for this system are readily ob- 
tained by extending the single-helicopter/single-cable 
system formulation in appendix B. The generalized co- 
ordinates include cable velocities for every cable, and 
the results represent either elastic or inelastic cables. 

The EOMs in the above applications work have been 
systematically formulated in terms of the natural vec- 
tors and matrices of rigid-body mechanics. This de- 
parts from the earlier formulations which were given 
in terms of the scalar components of all vectors. One 
consequence is that repreated matrix and vector prod- 
ucts are visible in the equation summaries, and their 
repeated calculation can be eliminated in efficient cod- 
ing. Another is that programming in a language that 
admits operators for such products would result in 
simpler, briefer code with associated savings in the 
programming task. These have not been significant 
issues for single-helicopter simulations, but they are 
of greater consequence for the dual-lift and multilift 
systems. 

5. APPROXIMATE NONLINEAR MODELS 

Approximations that reduce the number of DOFs 
and thereby simplify the system are always of interest 
in control analysis and design. These include the as- 
sumption of inelastic cables considered in the previous 
sections. Approximations of the rigid bodies as point 
masses or as rigid rods have also appeared in the litera- 
ture (e.g., ref. 10), particularly for dual-lift systems. In 
references 6 and 7, the dual-lift system with spreader 
bar is approximated by using a point-mass load and 
a rigid-rod spreader bar. In references 17 and 18 all 
bodies, including the helicopters, are approximated as 
point masses. Modification of the nonlinear EOMs to 
include these degenerate-body approximations is con- 
sidered in this section. 

The nature of these approximations in representing 
real systems is of interest. The rigid rod can approx- 
imate elongated bodies with negligible inertia about 
one axis. The point-mass assumption implies negligi- 
ble inertia about all axes. This is never realistic in 
slung loads, but it does approximate possible practical 
situations in which the attitude dynamics of a rigid 
body do not affect the remaining DOFs (reduced sys- 
tem). First, if the suspension is attached at the c.g. 
of a body and the applied forces are negligibly depen- 
dent on attitude, then the reduced system motion is 
independent of its attitude dynamics and is governed 
by EOMs identical to those obtained assuming a point 



mass. Second, if a helicopter can be assumed to control 
the slung-load system solely through its applied force, 
then its attitude dynamics are dependent DOFs as re- 
quired to generate that force, and can be removed from 
the system DOFs to leave a reduced system in which 
the helicopter’s applied force is the control. Simulation 
and control analysis of the reduced system can then 
proceed without the details of the helicopter model 
being considered, and separately from implementation 
of the applied- force controller in the helicopter. A cor- 
responding separation of the EOMs is examined below 
and it is shown to result in equations for the load- 
suspension subsystem that are identical to those ob- 
tained when point- mass helicopters are assumed. This 
separation of the slung-load control and the related 
point-mass approximations in references 17 and 18 are 
of interest for control synthesis based on inversion of 
the nonlinear model. 

Degenerate Body Approximations 

Point-mass bodies- If any bodies in the system 
are approximated as point masses, then all formula- 
tions of the EOMs (eqs. (8), (9), (14), (16)-(20)) can 
be applied by removing the attitude coordinates, an- 
gular velocities, applied and suspension moments, and 
inertias of the point-mass bodies from the configura- 
tion vectors and matrices r, v , X, /a, /c, and D. 

Rigid rods- A rigid rod has a singular inertia ma- 
trix so that equations (9) and (16) for the uncon- 
strained system, in which D~ l occurs, are invalid. 
Repairs can be made by defining a reduced configu- 
ration velocity that contains only as many coordinates 
as there are DOFs of the unconstrained system and 
retracing the derivations of section 3, beginning with 
the Newton- Euler equations. 

A rigid rod is a line segment in space with only two 
attitude DOFs. Attach body axes — {ib,jb,kb} 
to the rod with ib along the rod and located in inertial 
space by its heading and pitch angles: 


The overbar ( ) is introduced to indicate objects from 
the three-DOF rigid-body attitude dynamics, such as 
Wbb, ab , which are reduced to represent rigid rods. 
The direction angles are taken as heading and pitch 
here, with the assumption that the rod is never ver- 
tical. The columns of Wb are. the axes,of rotation 
for pitch and heading. Also, ab = Wb b uibb, since 

Wb T b Wb b = I ■ 

Next, consider a system of n bodies containing rigid 
rods, where 71 are the indices indicating the rigid rods, 
and wr, or are lists of their angular velocities and 
rates; a>r = (ujjj , j e 72.), af — ( aj , j e 11). To avoid 
the difficulties in deriving the EOMs from equation (7) 
that result from the singular £>, define the reduced con- 
figuration velocity v to contain the Euler-angle rates 
of the rigid rods. The relation v(y) is 



where 

Wr = diag{Wjj, j c 11} 

Here, v has been partitioned into the angular velocities 
of the rigid rods a;r and the remaining configuration 
velocities ul, and v contains 6n — r coordinates cor- 
responding to the DOFs of the unconstrained system. 
Also, v — W T v since W T W = L 

Next, let u be the generalized velocity coordinates 
for the system and denote the relations t?(w), v(u) as 

v — A u and v = A u 

where, from equation (32), A, A are related by 

A = WA and A = W T A 

Since v = W Au, replacing A in equation (8) by W A 
gives the following convenient form for the EOMs: 


Jb - diag{ 0, Jb\ Jb f } 

T biN = ( 0 b ) £ 3 (^ 6 ) 

0 — sin 0 b 

1 0 

0 cos 8 b 
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The configuration forces fo, fc are the usual objects 
in 6n-dimensional space previously defined, with forces 
and moments listed in the same order as the rigid-body 
linear and angular velocities in v above (eq. (32)). 

For the constrained system with inelastic cables, 
equation (13) expressing d’Alembert’s principle still 

applies ( A T fc = A T W T fc — 0), and then equa- 
tion (33) can be rewritten in a form analogous to equa- 
tion (14): 


u = [A T DA]- 1 A r W T fo (34) 

For the system with elastic cables, A is a nonsingular 
square matrix so that equation (34) can be expanded 
in a form analogous to that of equation (9a): 

u = A~ 1 W 1 W T (fo + fc) (35) 


The coordinates u can be selected as in equation (16) 
to contain 6 n — r — c coordinates ul, representing the 
configuration motion with inelastic suspension, and 
c coordinates A, which define the motion induced by 
cable stretching. Equation (35) then becomes 


u — 




D 1 W T (fo + fc) (36) 


For inelastic cables, the condition A = 0 yields c scalar 
equations in the constraint force. This force may be 
expressed as 

Yc=W T fc = Hs 

where s is c independent parameters and H is a basis 
of the 6n — r dimensional constraint force space for 
/c. As in section 3, it can be shown that A and H 
are bases of the same space, and then the equations 
obtained above yield the solution for s as 

s = -[H T D~ l H }- 1 H T D~ l W T fo (37) 

where H = W T H, and where H is some known basis 
in 6n-dimensional space, or A from the kinematics u(F) 
in (6n — r)-dimensional space. Equations (36) and (37) 
are analogous to equations (16) and (19), respectively. 


Applications- Example applications are given in 
appendixes D and E. In appendix D, the dual-lift sys- 
tem with spreader bar is approximated as three point 
masses, representing the two helicopters and load, and 
a rigid-rod spreader bar. This system has 14 and 
10 DOFs for elastic and inelastic cables, respectively. 
The result completes the work begun in reference 17 
to obtain a degenerate-body approximation for this 


system. In appendix E, equations are given for the 
multilift system with m helicopters and a pendant sus- 
pension with all bodies approximated as point masses. 
These extend the equations given In reference 18 for a 
two-helicopter, two-dimensional point-mass system. 


Center-of- Gravity Attachments 


To examine coupling between the attitude DOFs of 
body Bn in an n-body system with the remaining 
DOFs, partition the configuration vectors and matri- 
ces to separate its attitude coordinates, inertia, and 
moments: 


r — 




fa = 


fal 

MAn-n 


[ fd 

\MCn n 


The Newton-Euler equations (9b) for the configuration 
with elastic suspension are then similarly partitioned 
as 


i)\ = Dl 1 (fol + /cl) 

= Dl- l (fgl+fal-Xl+fcl) 
un n = Jn~ l ( MAn n — S(um n )Jn um n -f MCn n ) 

Coupling of the reduced system (rl,t;l,i;l) with the 
attitude dynamics (an, ion n , ujn n ) can occur in the 
aerodynamics, fal, and the suspension forces, /cl. 
Aerodynamic coupling is principally coupling of the 
body’s attitude and c.g. dynamics, plus secondary in- 
terbody interferences. This coupling can be neglected 



if the aerodynamic force is negligible (e.g., load aero- 
dynamics at hover and low speeds) or nearly indepen- 
dent of attitude (e.g., load aerodynamics when these 
are dominated by attitude-independent drag, or when 
the attitude is assumed to be constant). 

The bodies interact principally through fc. For an 
elastic suspension, fc can be given by equations (10) 
or (11) from the inertial cable directions and lengths, 
and the attitude angles of all bodies, and these are 
functions of the configuration position r. It can be 
shown that 

If Rn*j = 0, je Cn 

then MCn n = 0, /cl = /cl(rl) 

where Cn = indices of the cables attached to Bn at 
attachment points {Rn*j, j e Cn} on Bn 

Thus, if the suspension is attached at the c.g. of 
5n, and the dependence of /ol on its attitude is ne- 
glected, then the subsystem (rl y vl,vl) is independent 
of (an,ion n ,Con n ). Further, generalized coordinates 
for the reduced system can be selected according to 
the procedure of section 4 to represent both elastic and 
inelastic suspensions. The coordinates (rl, v\) are also 
those of a configuration in which Bn is assumed to be 
a point mass, and the Newton-Euler equations for this 
reduced system are identical to those obtained when 
Bn is assumed to be point mass. 

If Bn is a helicopter, then its aerodynamic force is 
attitude dependent. However, system control is sim- 
plified by a c.g. attachment at the helicopter. The 
helicopter controls vl solely through its applied force, 
jF0n;v, and attitude is selected to generate the required 
control force. For c.g. attachment, the helicopter’s at- 
titude dynamics are those of a free rigid body and 
are unaffected by the suspension or the motion of any 
other body. The extension of this separation of the 
control to more general single-point attachments is 
considered next. 

Simplified Helicopter Model 

Helicopter aerodynamics are a significant source of 
complexity in the simulation of multilift systems and of 
their control analysis and design. It is desirable to sep- 
arate the design of the helicopter control from that of 
the slung-load system control to allow a similar separa- 
tion of the simulation and analysis work. The general 
approach is to view the helicopters as force genera- 
tors that control the slung-load system by means of 


their controllable applied forces, and then to consider 
the slung-load control independently of the details of 
controlling the helicopters’ applied forces. The heli- 
copters’ applied forces are controlled through the heli- 
copter attitude DOFs; that is, given the instantaneous 
helicopter state and the desired force, the helicopter’s 
force and moment balance equations can be solved for 
the corresponding helicopter attitude angles and con- 
trol settings, and then the controls can be used to drive 
to its desired value. If the attitude control bandwidth 
is much faster than the bandwidth of the applied forces 
required to control the slung-load motion, then it can 
be assumed that the attitude is always approximately 
at the required value; that is, the attitude control is 
instantaneous and the attitude DOFs are in steady 
state (a; = 0) relative to the applied force variations 
needed to control the remaining DOFs of the system. 
The condition Co — 0 characterizes the “controlled he- 
licopter approximation” considered next. 

Partitioned EOMs- Let H be the indices identi- 
fying the helicopter rigid bodies { Bj , j e H} in an 
n-body slung-load system. Denote the helicopters’ an- 
gle coordinates, inertias, moments, and forces as a/i, 
a ;ft, Jh, /c2, /o2, foh , as defined in part a of ta- 
ble 2. Assume that 6n generalized velocity coordinates 
u have been selected for use with equations (16)-(20); 
that these include the helicopter angular velocities coh ; 
and that the remaining coordinates u 1 suffice to define 
the velocity of any point in the system outside the he- 
licopters. A review of the generalized coordinates se- 
lected for the examples in the appendixes show's that 
these are “reference-point coordinates”; that is, they 
contain the velocity of a reference point and additional 
coordinates that define the motion of all points rela- 
tive to the reference point. As such, coh is present in u 
for all examples, but the use of a helicopter c.g. as the 
reference point in most examples violates the above 
requirement for ul. This can be repaired by moving 
the reference point to the load, as in appendix E, or to 
a suspension attachment point on a helicopter. Next, 
partition the kinematics, the mass-inertia matrix, and 
the configuration forces u, v, v(u), u(u), D, /o, /c, 
as shown in part b of table 2, to separate the EOMs 
into equations for Coh and for the remaining (reduced 
system) DOFs, ul: 

ul = [Ain AI^D- 1 (fo + fc) 

= Ain Dl -1 (/ol + fcl)t + AI 12 Coh 

Coh — JTi" 1 (/o2 + /c2) 

In general, ill is coupled with the helicopter attitude 
dynamics (ah, coh, Coh) through the term A 1 12 Coh 
and the possible dependence of /cl, /ol, and Ain on 
(ah, coh). 
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Table 2. Equations of motion for controlled helicopter approximation 

(a) Helicopter attitude variables, inertias, moments, and forces: 

H = indices of helicopters in an n-body system: 

ah = ( aj , j e H) 

wh = ( u>jj , j € H) 

Jh = diag{Jj, j e Ti) 
fc2 = (MCjj, j e H) 
fo2 = (MOjj, j e H) 
foh = (FOjjVi j f hi) 


(b) Partitioned equations: 



For inelastic suspension: 

fc 


/cl 

fc2 


A1 

A2 


s = -[Al r D l" 1 A1 + A2 r Jh- 1 A2]- 1 (A1 T £>1 -1 /ol + A2 r Jh- 1 fo2) 


ii\ = AI U Dl- 1 (/ol + /cl) + AI 12 uih 
Coh = Jh~ x (/o2 + fc 2) 
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r; 


Table 2. Concluded 


(c) Controlled helicopter approximation: 

If: wh = 0 

Then: ill — A [ / Dl~ l (fol + /cl) 
fol = — /c2 

5 = -[A1 t Dl~ l At]” 1 A1 T Dl~ l fol 

If also: /cl = /cl(rl), An = An(rl), /ol = /ol(rl, ul , fob) 
where: rl = rl (rl, ul) and foh are the helicopter applied forces 
Then: (1) The reduced system EOMs are independent of the helicopter 
attitude DOFs: ill = ill (rl, wl, foh ) 

And: (2) The steady-state controls and attitude of conventional helicopters 

satisfy the following quasi-static trim equations: 

For j e H : 

Given: F0j N , aj, Vj* N , Jpj 
Find: (Sj, 4>j, Oj ) such that 

V A j* j = T jiN (Vj* N - W0 N ) 
ujj = Wjjiaj) 

T N j FAjj u)jj, Sj) + mj g N = F0j N + mj Aj u 

MAjj(VAj *j, ujj, 6j ) = S(u)jj) Jj wjj - MCjj 

where 

FOJn is an element from foh 

Aj are rows of A which give Vj* N = Aj u 

MCjj is an element of /c2 = A2(r) s(rl) 

WOn is the mean wind velocity 

Wjj defines angular velocity from Euler-angle rates of Bj 
(appendix A, eq. (72)) 
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Controlled helicopter approximation- If the 

condition ojh = 0 is imposed, then the coupling of 
ul with Coh drops out, and the helicopter applied- 
moments satisfy fo2 = —fc2 = — A2 s. Substitute this 
expression for fo2 in the equation for the interaction 
force parameters s to obtain the result shown in part c 
of table 2: s = -[Al T Dl~ l Al] -1 Al r Dl -1 fol. 
These parameters are independent of the helicopter 
moments and inertias, /o2, Jh. In general, A depends 
only on position variables. If Al depends only on 
position variables rl, which can be determined from 
ul, then /cl = Al s is independent of the helicopter 
attitudes, except for possible dependence of fol on 
(ah, ujh). This condition is met in all the applica- 
tions examples presented in the appendixes for the 
coordinates u selected there. In these examples, Al 
depends on the cable directions and the attitude an- 
gles of the bodies other than the helicopters, and ul 
contains the cable velocities and the angular velocities 
of these bodies. 

To examine /ol, partition it into the helicopter ap- 
plied forces and inertia reactions, foh ~ FOJn , jt H, 
and the forces and moments on the remaining bodies 
of the system (load, spreader bar,. . .): 

(F0j N , j i H\ 

f °Y = 

\M0jj, jin) 

( mj g N + FAj N - mj Aju j <£ H \ 
MAjj - s{vjj) jj ujj jin) 

The helicopter forces foh are considered to be the con- 
trols of the reduced system. The remaining forces, 
/ol 7 , contain the aerodynamics and the inertia reac- 
tions of the remaining bodies. Coupling of these aero- 
dynamics with (ah, ojh) can be neglected. In the iner- 
tia reactions, Aj are the rows of A that define V j* N (u) 
and that can introduce coupling with ( ah , a )h). In 
the examples in the appendixes this occurs as centrifu- 
gal accelerations in the helicopter attachment- point 
moment arms, but these terms are zero for straight- 
line flight (u)h « 0) and are expected to be negligible 
otherwise. 

Last, the submatrix AIu = [V^ ul] is independent 
of ah in the examples in the appendixes for the coor- 
dinates u selected there. 


Reduced system EOMs- For the examples 
in the appendixes, the reduced slung-load system 
(rl, ul, ul) can be approximated as being inde- 
pendent of the helicopter attitude variables by using 
the controlled helicopter condition and neglecting sec- 
ondary coupling in foV. The helicopter applied forces 
foh are the controls for the reduced system. 

Examples are omitted since the reduced system 
equations are identical to the equations obtained by 
assuming point-mass helicopters. For point-mass heli- 
copters, (jjh, /c2, and fo2 are undefined, and for c.g. 
attachments at the helicopters, AI X2 = 0 and fc2 = 0, 
so that the resulting equations for ul and s are identi- 
cal, in both cases, to those of the reduced system given 
above. 

Quasi-static equations for helicopter attitude 
and controls- If the reduced system is known, then 
the helicopter attitude angles and controls are depen- 
dent variables that can be determined by solving the 
quasi-static trim equations given in part c of table 2. 
The helicopters are assumed to be conventional, with 
four controls, and F0j Ni Vj* N , aj are known for each 
helicopter, Bj , jeTi , along with a specified value of 
the redundant attitude angle tpj. Then the four heli- 
copter controls and remaining two attitude angles can 
be determined from the kinematic equations and the 
force- and-moment balance equations listed in part c 
of table 2. The quasi-static trim equations are related 
to the usual static trim algorithm found in conven- 
tional helicopter simulations. The helicopter model 
(FAjj, MAjj ) has the usual independent variables as 
noted in the table. For isolated helicopters in static 
equilibrium, the angular velocity, linear acceleration, 
and suspension forces and moments are zero. Here, the 
equations account (1) for nonzero angular motion, aj , 
which is presumed to be known from the attitude his- 
tory; (2) for linear acceleration and suspension forces 
in F0//vs which are known from the reduced system; 
and (3) for the suspension moments, MCjj , in the 
moment equation. The body axis components of the 
suspension moments and rigid-body velocities depend 
on attitude, and these relations are included in the 
iterative solution procedure. 

6. LINEARIZED EQUATIONS OF MOTION 
FOR GENERAL SLUNG-LOAD SYSTEMS 
WITH INELASTIC SUSPENSIONS 

Linear analysis is the most commonly used technique 
in stability studies and control system design, hence 
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linearized equations of motion (LEOMs) for slung-load 
systems are of interest. These define the perturbation 
behavior in the neighborhood of a single flight condi- 
tion, and this behavior is characterized by the eigen- 
vectors and eigenvalues of the linear equations. 

The LEOMs can be obtained by numerical evalua- 
tion of the gradients from a working nonlinear simu- 
lation program (e.g., ref. 15) or by analytical deriva- 
tion from nonlinear EOMs. Previous derivations in the 
slung-load literature are specific to the suspension and 
flight conditions treated, and to the approximations 
made. Early studies of single-helicopter systems con- 
sider single and multipoint suspensions with various 
simplifications of the load and its aerodynamics (refs. 9 
and 10), or without such simplifications (ref. 33). More 
recent studies have focused on transport of the MIL- 
VAN and bluff bodies whose stability as slung loads 
is problematical (refs. 11, 15, and 16). The LEOMs 
are derived from both Lagrange’s equations and from 
the Newton-Euler equations, and inelastic cables are 
assumed in this literature. 

Linearized equations of motion for the twin-lift sys- 
tem with spreader bar are derived from Lagrange’s 
equations in references 6 and 7 for hover conditions 
and point-mass load without load aerodynamics. Ad- 
ditional linear control studies reported in references 37 
and 38 utilize these results. 

The present work provides formulations of the 
LEOMs for general slung-load systems at general ac- 
celerating or static equilibrium reference flight con- 
ditions. The cables are assumed inelastic, and con- 
ventional aerodynamics for all bodies are included. 
Controllable suspension parameters are omitted for 
brevity, but the derivation of the LEOMs can be ex- 
panded routinely to include these. The LEOMs are 
derived from Lagrange’s equations (eqs. (23) and (24)) 
and from the Newton-Euler formulations of the nonlin- 
ear EOMs (eqs. (14), (19) and (20)) previously given. 
In appendix F, the general formulation from equa- 
tion (14) is applied to the dual-lift system with a 
spreader bar. Results are given for static equilibrium, 


and account for load and spreader-bar inertias and 
aerodynamics, and for dissimilar’ helicopters. 

The derivation in each case is carried out in three 
steps. First, a general form of the nonlinear equations 
showing the dependencies of its terms on the indepen- 
dent variables of the problem is stated. These depen- 
dencies represent typical working simulation models. 
Second, the reference trajectory or flight conditions 
are defined. Third, the perturbation equations are de- 
rived and evaluated on the reference trajectory. In all 
cases, a representative model of the applied forces and 
their linearization and a definition of reference trajec- 
tories are needed; these matters are settled first. 

Applied External Forces 

The applied external forces and moments on each 
body Bi are due to gravity and aerodynamics: 

Fix = mi gN + T^ ti F Aii ( V Ai **, on*, Si) 

Mii = MAii {VAi\, oh* , Si) (38) 

where 

VAi\ = T un {Vi* N - W0 N ) 

It is assumed that the aerodynamic and rotor forces 
and moments of Bi are given in simulation models in 
body axes as functions of the body axis coordinates 
of the c.g., the velocity relative to the air mass V Ai* i 
(WO is the mean wind velocity), the angular velocity, 
uii , and controls, Si. This general, conventional model 
neglects position- and acceleration-dependent aerody- 
namics (e.g., ground effects, other altitude- dependent 
effects, interbody interference effects, unsteady aero- 
dynamics), and gust disturbance effects, but is ade- 
quate for most loads, systems, and flight conditions. 
Most of the omitted effects can be superposed or 
are only weakly position-dependent. Models for the 
higher-order unsteady cargo carrier aerodynamics are 
discussed in references 11, 15, and 16. 

Equation 38 can be assembled for the complete con- 
figuration as 
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(396) 


The notation T, fba , w, va is conveniently introduced 
above, where /6a, va contain the body axis compo- 
nents of the aerodynamic forces and c.g. velocities. 
The LEOMs require the gradients of the aerodynamic 
and rotor forces: 

Fbva = v£, fba 

Fb s = Vj fba (40a) 

where 5 is a list of all controls (helicopter rotor con- 
trols and control-surface deflections). The terms Fb va 
and Fb s are assembled from the stability and control 
derivatives for each body: 

( FAU \ 

(^vaui Vji) j i = l r ..,n (406) 

\MAiJ 


accelerating straight-line flight, and various kinds of 
turns) connected by brief transitional maneuvers. For 
multibody slung-load systems the relative motion of 
the bodies is usually carefully controlled such that the 
configuration moves as a single rigid body during each 
steady segment, or very nearly so. The relative posi- 
tions of the bodies differ from one segment to another, 
so that the configuration must be rearranged during 
the transitions between these segments. In the follow- 
ing discussion, the term “fixed configuration” means 
that the system behaves as a rigid body. That is, if 
Rl* is a reference point in the configuration, and Rj is 
any other point in the configuration, then the inertial 
velocity of Rj is 

Vj = Vl* ~f u; x Rl*j 


These derivatives are all assumed to be available from 
a simulation model for any flight condition of inter- 
est. Nonlinear and linear models for helicopter aerody- 
namics are commonly available, but the aerodynamics 
of loads are much less systematically known or pre- 
dictable. Static aerodynamic data for several military 
vehicles and cargo carriers are reported in literature 
cited in reference 39. Stability derivatives are less well 
known, but reference 33 contains data for a vehicle, 
and references 11, 15, and 16 review the literature and 
data for the MILVAN and cargo carriers. 


where uj is the configuration angular velocity. If a; = 0 
then the configuration is nonrotating as well, and every 
point has identical velocity. This occurs in hover, in 
static equilibrium, and, approximately, in accelerating 
straight-line flight. 

Linear analysis is usually limited to hover and static- 
equilibrium flight conditions. Static equilibrium is de- 
fined here as a fixed, nonrotating configuration that 
follows an unaccelerated straight-line path: 

v 0 (t) = 0 


Reference Trajectories 


f*\ 


The linearization is taken about a reference configu- 
ration trajectory (r 0 (£), v 0 {t), v Q (t), 6 0 (t ), t 0 <t<t f } 
which satisfies the EOMs and additional engineering 
conditions which specify redundant variables. Quanti- 
ties evaluated on the reference trajectory are indicated 
by the subscript ( ) G . 


Vo 


I 

0 


V0 N 


Vo/ 

— v 0 {to) ~b v 0 (£ £ 0 ) 


(41) 


For a single rigid-body aircraft, a complete reference 
trajectory consists of a sequence of extended quasi- 
static flight segments (hovering, static equilibrium, or 


6 0 (t) = 5 0 (t 0 ) ~ constant 

where VO is the fixed reference velocity. Hover is the 
special case of equation (41) with no motion, VO = 0. 
The quantities (r 0 (£ 0 ), <5 0 (£ 0 )} are obtained from a 


26 



trim-solution algorithm. References 15 and 16 describe 
trim algorithms for single helicopters with one- and 
two-point suspensions. No complete algorithm has 
been given yet for any dual-lift system, but a study 
of equilibrium configurations along general reference 
trajectories for the dual-lift system with spreader bar 
obtained by solving a simplified trim problem is de- 
scribed in reference 5. Results are given for configura- 
tion geometry, cable forces, and thrust requirements. 

The reference trajectory generalized coordinates 
Qo{t), u 0 (t),u 0 (t) are needed for the linearization, and 
are obtained from the reference configuration motion 
using the following geometry and kinematics: 

r = r(q) 

v — A(q) u + w (42) 

v = A u + A u 

Here, the kinematic model assumes that the coordi- 
nates u are chosen to generate configuration velocities 
relative to the mean wind, v — w, as is usual for lin- 
ear analysis. If v = 0, then Au = — Ail, and it does 
not follow that u is fixed in static equilibrium with- 
out further assumptions about the generalized coor- 
dinates. If u is selected such that A(q) depends only 
on variables that are fixed in static equilibrium, then 
A{q o ),u 0 are fixed and A 0 ,ii 0 are zero in static equi- 
librium. Such variables can be the coordinates of any 
line segment connecting points in the configuration, 
such as the rigid-body attitude angles and cable an- 
gles relative to inertial or body axes. The coordinates 
selected for the nonlinear EOMs in all the applications 
cases presented in appendixes B-E are reference point 
coordinates of this type. That is, (#, u ) contain the po- 
sition and velocity of a reference point V1* N ) 

and additional coordinates (g,rl), which suffice to de- 
fine the positions and velocities of all points in the 
system relative to the position and velocity of the ref- 
erence point 

(m* N \ (vv N \ 

<7=1 I , u = I I , v = A(q) u 

where q is constant and u = 0 if the configuration 
is fixed and nonrotating. For linear analysis, the 
reference- point coordinates are usually given relative 
to the mean wind in body axes: 

/ JVAll dt\ /VAVA 

<7=1 J , u = ( 1 , u = A(q) u + w 


wdiere, as above, q is constant and u = 0 for fixed, non- 
rotating configurations. Then, for static equilibrium, 

Qo(t) = q(r 0 (t)) 



Linearized Equations from Lagrange’s 
Equation 

The nonlinear EOMs are obtained from Lagrange’s 
equations in terms of the generalized position coor- 
dinates q . Equations for general slung-load systems 
are given by equations (23) and (24) and these are re- 
peated in table 3. The terms Af, k are obtained from 
derivatives of the kinetic energy. In the generalized 
forces, A defines v{q,q), and T,fba,va are defined in 
equation (39). All dependencies on the independent 
variables {q, q, q,8] are indicated. 

The reference trajectory {q 0 {t) y q Q (t), <7ofa)? <5o(£)} 
satisfies the EOMs. The LEOMs are obtained by ex- 
panding these to first order in variations about the 
reference trajectory: 

q = q 0 + 8q 

Q = Qo + 8q 

Q = Qo + 8q 

8 — 8 0 + A 

and then subtracting the EOMs evaluated on the ref- 
erence trajectory to get 

M 0 8q + 6M q Q + [V^] 0 5q + [Vp;] 0 8q = 8Q 

The variations 8 A/, 8Q are further expanded by using 
their general dependencies on the independent vari- 
ables of the problem. Complete results are given in 
the equation summary, table 3. Products and vari- 
ations are interchanged where useful, for example, 
8M(q) q 0 = 8{M{q) q Q ) = [Vj (M (q) q 0 )} qo 6q. 

The terms in the perturbation equations that are 
due to the dynamics M (q) q and k(q,q) are Af OJ A/ 0 , 
and k q , k q . The term M 0 is zero whenever q Q — 0; 
the coordinates q are usually selected such that this 
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Table 3. Linearized EOMs for slung-load systems from Lagrange’s equations 


Nonlinear EOMs: 

M ( q ) q + k(q, q) = Q(q, q, 6) 
v = A(q ) q + w, w = 0 
KE=\ v T D v 
M = Vg KE = A T D A 

k = [Vj V 4 KE\q - V q KE = A T DA q + (A - G) T D v 
G = [V^ (A(q) q)} 

Q = A(q T ) ( fg + fa) 
fa = T(q) fba(va , S) 
va = T(q) T ( v — w) 

LEOMs: 

Reference trajectory: {q 0 (t), q 0 (t), <io(t), S Q (t)} 

Perturbations: {6q, 6q, 6q , A} 

Equations: M a 6q + C Q 6q + K 0 6q = Qs A 

where C a = k q — Q q 

K 0 = Mq + kq — Qq 

Qs = Af T 0 Fbs 

and Qq = Af F vaN A 0 M 0 = [V^(M(q) q 0 )\ Q , t 

Qq = Q^-q + Q^q kq = ^'qk\q 0 ,q 0 

Ql q = [VT(A(q) T (fg + T(q) fba 0 ))] qo k q = \V T q k) qnAo 
Q2 q = A T 0 T a V q fba = A T 0 F vaN [V^ (T a T(q) T (v a - tit) + A(q)q 0 )} ^ 
where F vaN — T a Fb va Tj ; Fbs, F va , T, fba are defined in eqs. (39) and (40) 
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occurs in static equilibrium. All terms in k q , k q contain 
a velocity coordinate from q 0 and are zero whenever 
Qo = 0; the coordinates q are usually selected such 
that this occurs in hover, and then k(q,q) need not be 
formed, and 

6q = M 0 *( Qq 5q + Q q Sq + Qs A) (44) 

The terms that are due to k also drop out or simplify in 
static equilibrium for special choices of the generalized 
coordinates as established in the next subsection. 


A reference trajectory is any {$,(*), u Q {t), u 0 (t ), 
S 0 (t), to < t < tf} which satisfies the EOMs. The 
LEOMs are obtained by expanding the EOMs to first 
order in the independent variations. { 8q , 8u , 8ii , A}, 
and then subtracting the EOMs evaluated on the ref- 
erence trajectory: 

M 0 Six + 6M u 0 + 6A T (X 0 + D£ 0 ) 

+A^{SX + DSC) = SQ 


The terms in the LEOMs that axe due to the applied 
forces are Qs, Qq, Q q . To obtain these expressions, the 
extended chain of dependencies 

= Q(A(q),fa(T(q),fba(6 } va(v(q,q)))) 

should be noted. It can be seen in table 3 that Qs con- 
tains the control derivatives and that Q q contains the 
stability derivatives. Furthermore, Q q has been sepa- 
rated into Ql q , which contains the equilibrium forces, 
and Q2 q , which contains a combination of stability 
derivatives with the reference velocity. The term Q2 q 
is zero in windless hover where q 0 ,v 0 — w are zero. The 
gradients of T, A in Q q expand principally to gradi- 
ents of coordinate transformations in the applications 
work; useful general formulas for this are given in ap- 
pendix A (table 7). 

Linearized Equations from d’Alembert’s 
Principle 

The nonlinear equation (14) can be arranged in a 
form that parallels the result from Lagrange’s equa- 
tions, and this form is given in the equation summary, 
table 4. All dependencies on the independent variables 
{q,u,ii,8} are noted. These EOMs differ from the La- 
grangian equations in the formulation of the second- 
order velocity terms k(q,u) and in the use of gener- 
alized velocity coordinates, u different from q. The 
second-order velocity term X is formulated in equa- 
tion (3) in terms of the system’s angular velocities, 
which are coordinates of t>, whereas A is formulated 
in the applications work (appendixes B-E) in terms of 
{• Q » u } in all cases; these dependencies are assumed in 
table 4. The notation £ = Au, introduced for conve- 
nience, represents a configuration acceleration. Equa- 
tion (14) is applied by formally defining the general- 
ized velocity coordinates u, and position variables are 
introduced as needed to define v{q,u) and T(q). In 
general, fewer than d such variables may suffice for 
this purpose, and these are assumed to be given by 
a linear relation q = Y(q) u in table 4, where Y has 
full rank. If q is expanded to a set of d coordinates 
sufficient to define r(£), then Y is nonsingular. 


The variations 8M, 6 A, 6X, <5f, SQ can be expanded 
further; the results are given in the equation summar 3 r , 
table 4. In the coefficient matrices of the linearized 
equations, M 0) M 0 are from the term M(q) u ; k u ,k q 
are from k(q, u ); and Q u , Q q , Qs are from the general- 
ized forces. 


In static equilibrium, the rigid-body angular ve- 
locities are zero, so that X Q ,X v ,X u ,X q are all zero 
regardless of the choice of generalized coordinates. 
Furthermore, the generalized coordinates used in all 
application examples given in appendixes B-E are 
reference-point coordinates. If inertial coordinates of 
the reference-point velocity, FI* N , are used as in these 
appendixes, then it is readily shown that the terms k u 
and k q due to k , and M 0) are zero in static equilib- 
rium. However, for linear analysis it is customary to 
use body axis coordinates of the velocity relative to the 
mean wind, V Al* i . In this case, A has the partitioned 
form 
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where q is fixed and u = 0 is in static equilibrium, and 
where A depends only on coordinates (q, cel), which 
are fixed in static equilibrium. Then 


M 0 = 0 
&> = 0 
^ = 0 


£« = Vj Au 0 - -V£ 
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Table 4. Linearized slung-load EOMs from equation (14) 


Nonlinear EOMs: 

M(q ) u + k(q,u) — Q(q,u,S) 
q = Y(q) u 

v = A(q) u + w, w = 0 
M = A(q) T D A(q) 
k = A(q) T {X(v) + DZ(q t u)) 

S = Mq’ u ) u 

Q = A(q) T (fg + fa) 
fa = T(q) fba(va y 5 ) 
va = T(q) T ( v — w) 

LEOMs: 

Reference trajectory: (g 0 (t), u G (f), u 0 (£), 5 0 (t)} 
Perturbations: <5u, A} 

Equations: M 0 5u + C Q &u + 6q = A 
Sq — Y 0 5q — Y 0 Su 

where 

F 0 = [V^ Y(q) u 0 ] qo 
O o k u Qu 

K 0 = Af 0 + k q — Qg 

Qi = T o Fb 6 

and 

Mo = [V^( M{q ) n 0 )], o 

= iC I A 'u + D e„] 

A:, = A t 0 [X q + D e,] + [V^4(g) 7 ’(X 0 + D 
where 


= [vprk 
X u = Ao 
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Table 4. Concluded. 

X 9 = JffA(q) u 0 } q „ 

= Ao + [VlA(q 0 ,u) Uo] 

6, = [V^ A(q, u 0 ) u 0 ] qo 
and 

Qu ~ F yaN 

Qq = Q^-q + Q%q 

Ql g =M{A(q) T (f9 + T(q) fba 0 ))] Vo 
Q2 q = A T 0 F vaN [Vj (T 0 T{q) T ( v a -w)+ A{q) u 0 )} Qo 
where F vaN = T a Fb va Tj ; Fb s , F va , T, fba are defined in eqs. (39) and (40) 
Second-order ODE for 8q: 

M' 0 6q + C' 0 6q + K' 0 6q = Q s A 
u = U{q) q 
K = M 0 U 0 

C' Q = C 0 U 0 + Mo {Uo + (U(q 0 . q) q 0 )) 

K' 0 = K a + C 0 [y^U(q)q 0 ] g<1 + M a [Vj {U ( q , q 0 )q 0 + U(q) Qo)] qo 
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where VAO = VO — WO. Assuming that u>li is in- 
cluded in the generalized coordinates u, then the only 
nonzero column of is the gradient with respect to 
wli: 


/A 


6u = Sufl\ 


I 

0 


S(VA0i)6u>li 


Vo/ 


(46a) 


For this case, the LEOMs for static equilibrium sim- 
plify to 

5ii = M 0 (( Q u — D£ u )6u + QqSq + Qs A) (466) 


Here, £ u is the only contribution from the second-order 
velocity terms k. It is given from above; it is zero in 
hover, and otherwise adds terms in VAO x 6cul to the 
velocity coefficient matrix in the LEOMs* The equa- 
tions for static equilibrium and hover differ in form 
only by the terms and Q2 q from Q q \ these terms 
are proportional to the reference airspeed. 

The LEOMs can also be given as a second-order ordi- 
nary differential equation (ODE) in the generalized po- 
sition coordinates. This equation is obtained by form- 
ing the variations Su(6q,6q), 6u(6q,6q,5q ) from the 
kinematic relation u — U (q) q and substituting these 
in the state equation. Here, q contains d coordinates, 
which suffice to define r(£), and U is nonsingular. Re- 
sults for general reference conditions and coordinates 
are included in table 4. For static equilibrium and 
reference-point coordinates such that 


Linearized Equations from Explicit Constraint 
Method 

The nonlinear equations (19) and (20) for the in- 
elastic system are repeated in the equation summary, 
table 5. This form assumes that 6n coordinates u 
have been selected for the elastic system, which con- 
tain 6n-c coordinates, ul, of the system with inelastic 
suspension, and c coordinates, A, which suffice to de- 
fine the motion resulting from cable stretching. The 
matrices A, A -1 define the relations v(q,u ), u{q,v), 
respectively, and the matrices Al, All , A, which ap- 
pear in the EOMs for the inelastic subsystem ul, are 
submatrices of A, A -1 . Equations (19) and (20) dif- 
fer from equation (14) in the formulation of the co- 
efficient matrix of fo, which is known analytically in 
equations (19) and (20) up to the inverse of the c x c 
matrix, S. The EOMs are linearized by expanding to 
first order in the variations {Sq, Su, Sii, A} and then 
subtracting equations (19) and (20) evaluated on the 
reference trajectory {q 0 (t),u 0 (t),u 0 (t),6 0 (t)}: 

6iil=SAIl T Al 0 ul 0 + AIlZ D~ l [Sfo + Sfc] (4 8a) 

In the first term, D~ 1 (fo 0 -f fc Q ) has been replaced 
with its equivalent, Al 0 iil 0 . The variation of fo is as- 
sembled from the variations of its terms, all of which 
were previously treated in tables 3 and 4. The con- 
straint force perturbation is obtained from variations 
of 

fc = E fo = A s 
s = — S _1 A t D~ l fo 
S = A t D~ l A 

which yield 

Ss = -S-'ISS s a + 6A t D - l fo 0 + AjD- 1 Sfo ] 


u = 


(TO 



Vo 


0 ^ 

U(q)J 



where 


SS s 0 = 8A t D~ l f Co + Aj D~ l SA s Q 
and then 


where u, q have the same meaning as above and are 
zero in static equilibrium, the second-order ODE form 
simplifies to 

M 0 U 0 8(j-(Q u -Ao D £ U )U 0 Sq-Q q Sq = Q s A (47) 

In applications, U is usually block diagonal and 
the modifications of the coefficient matrices in equa- 
tion (46b) that are required to obtain the coefficient 
matrices in equation (47) are simple. 


Ss = -S-'tf A T A 1„ «1 0 + AJ D- 1 8A s a Aj D~ l Sfo] 

(486) 

and also 

Sfc = <5 A s 0 + A 0 Ss = E 0 6 fo 

+ (7 4- E 0 )6A s 0 — A 0 S^ 1 6A t Al 0 ul Q (48c) 

In equation (48b), AI 0 ul Q replaces its equivalent as 
above. The completed expansions of the constraint- 
force perturbations 6s, Sfc are included in table 5. Fi- 
nally, the results for Sul in table 5 can be obtained by 
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Table 5. Linearized slung-load EOMs from equations (19) and (20) 


Nonlinear EOMs: 

v = A(q) u + w = Al ul + L A -f- w, 

ful\ ( AIl T 

I j = A~ l {q) (v-w) = f ^ 


w = 0 


(v — w) 


For A - 0 : 

Q = Y{q) 

(/o + /c) 

/c = E(g) /o = A(g) s 
E = -A(g) S(q)~ l A(q) T D~ l 
s = -S(q)~ l A(q) D~ l fo 
S = A(q) T D- 1 A (q) 

fo = fg + T(<?) /ba(va,£) - X(u) - jD £(<?,ul) 
ua = T(g) r (u — w) 

£ — Al(q , ul) ul 

Linearized suspension forces (A — 0,u; = 0): 

6s — <5^ H- <5ul “I - 5^ A 

Sfc = A a 6s+ [V^ A (q) s 0 ] qo 6q 
S q = Sq + 5 0 [i\j — D — Xq\ 

S u — So[FvaN Al 0 D Xn\ 

s& = So To Fb & 
and 

So = -So 1 A T o D 

S q = A {q) T A\o ul 0 ]q 0 
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Table 5. Concluded. 


LEOMs: 


6ui = r q 6q + r u 6ui + r s a 

Sq = Y 0 Sq -f- Y 0 Sul 
where 

Yo = [VJ Y{q) ul 0 ] qo 

Fq~ T q Go[Fq — D £ q — X q \ 

r u = Go[F vaN Alo — D £ u — x u ] 

F$ — G 0 T 0 Fb$ 

and 

G 0 = AIll D~ l [/ + E 0 ] 

X u = X v Al 0 

X q - [V^ Al(q) ul 0 ] Qo 

*v = [V? X} Vo 

£u = Al 0 d~ |V£ A\(q 0 ,vX) ul 0 ]ul 0 

£7 = Al{q,ul 0 ) ul 0 } qa 

and 

T q = V q [AIl{q)T - AIlJ D~ l K 0 S~ l k{q) T } qo Al a ul 0 
Fq = Fl q + F2 q 

Fl q = [V g r (T(q) fba 0 + A(q) So )]^ 

F2q = F va N [V^ (T a T(q) T (v 0 — w) + Al(q) nl 0 )] 9o 
where 

So = -S- 1 A T 0 D~' f 0o 

FvaN = To Fb va Tj 

Fbs , Fb va , T, fba are defined in equations (39) and (40) 
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combining equations (48a) and (48c) and expanding 
the remaining variations. The gradients of the second- 
order velocity dynamics X, £ and the stability and con- 
trol derivatives F va N,Fs occur in the coefficients of 
5u,6q, A as previously seen in table 4. The matrix 
r q is linear in the acceleration coordinates ul G , and is 
zero whenever ul Q — 0. The reference trajectory forces 
occur in F \ q , and differ from their occurrence in the 
analogous Ql q of table 4 in that a term in the cable- 
tension parameters s replaces a term in the external 
forces. 

Finally, if the coordinates, q t u, are reference- 
point coordinates with body axis coordinates of the 
reference-point velocity (eq. (45)) then for static equi- 
librium X 0 , X u , X q , r q are zero, is given by 
equation (46), and the LEOMs simplify to 

Sul = G 0 [{F vaN Al 0 - D £ u )6ul + F q 8q + T 0 Fb 6 A] 

(49) 

The matrices G a , F q are defined in table 5. 

Linearized Equations for Dual-Lift Systems 

Linearized equations of motion for the dual-lift sys- 
tem with spreader bar are derived in appendix F for 
static equilibrium by using the formulation from equa- 
tion (14). Reference- point coordinates are used with 
body axis components of the reference-point velocity, 
and equation (46) is applied. These results extend the 
hover equations given in references 6 and 7 to general 
static-equilibrium flight conditions, and include load 
and spreader-bar inertias and aerodynamics, and dis- 
similar helicopters. 

The coefficient matrices are expanded to a working 
form in terms of the natural vectors and matrices of 
three-dimensional rigid-body mechanics. The number 
of such terms is much greater for the linearized “equa- 
tions for 6u than for the nonlinear EOMs for u , but the 
derivation is feasible with the methods of this report. 
MACSYMA was used to expand matrix products in 
order to reduce error probabilities in the analytical 
expressions. 

7. CONCLUSIONS 

Nonlinear simulation equations for general slung- 
load systems have been derived. These account for 
any suspension geometry, including controllable ge- 
ometry and both elastic and inelastic suspensions, for 


any number of helicopters, and for any choice of gen- 
eralized velocity coordinates. Two formulations are 
given that generalize the previous case-specific con- 
ventional formulations for elastic and inelastic suspen- 
sions that are given in the slung-load literature. These 
formulations differ significantly in form and computa- 
tional requirements. A third, new formulation is given 
by selecting the generalized coordinates of the uncon- 
strained system to represent the constraints of the 
inelastic suspension. The internal suspension forces 
are calculated explicitly for both elastic and inelastic 
suspensions. The new formulation improves computa- 
tional efficiency significantly for inelastic suspensions, 
enables integration of elastic and inelastic suspension 
models in a single equation set, by using a single set of 
generalized coordinates, and is readily applied to the 
complex dual-lift and multilift systems. 

An inelastic suspension imposes constraints on the 
system motion, but their number is small compared 
with the number of DOFs. This fact, along with the 
cable interbody connections, distinguishes the slung- 
load systems from typical applications considered in 
the literature on multibody systems, and accounts for 
the efficiency of the new formulation for inelastic sus- 
pensions compared with the conventional ones in the 
slung-load and multibody literatures, the latter being 
efficient for highly constrained systems. 

In past work, simulations have usually utilized elas- 
tic suspension models and rigid-body velocity coordi- 
nates, whereas control analyses have been based on 
the inelastic suspension model. These methods are se- 
lected because of the analytical simplicity and compu- 
tational efficiency of the former, and the elimination of 
DOFs characterized by small motions in the latter. For 
simulation, the present results allow the use of any gen- 
eralized velocity coordinates with the elastic suspen- 
sion model, and the use of inelastic suspension mod- 
els with reduced computational penalties. For control, 
the formulation of efficient nonlinear EOMs for inelas- 
tic suspensions makes it possible to apply the recent 
global inverse-model methods to slung-load systems, 
and facilitates the derivation of linearized EOMs for 
single-flight-condition designs. 

Application of the general equations to the deriva- 
tion of simulation equations for specific systems is 
demonstrated for a series of single-, dual-, and 
multiple-helicopter systems. Results are given in pro- 
grammable form, with the dynamics formulated in 
terms of the natural vectors and matrices of three- 
dimensional rigid-body mechanics. This formulation 
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allows the devices of efficient coding to be applied to 
the vector-mechanical structure of the equations. 

The single- helicopter systems with single attach- 
ment points at the helicopter are readily treated, and 
the results provide alternative formulations of the 
EOMs to those already given in the literature. Re- 
sults for single-helicopter systems with multiple at- 
tachment points at the helicopter are more difficult 
to derive and were given solely for the bi filar suspen- 
sion, which has also been treated in the literature. 
Nonlinear equations for the inelastic inverted-V and 
inverted-Y suspensions are not found in the open lit- 
erature. The methods of this report do not address the 
analytical difficulties of these cases, which are omitted 
from the present work. 

The principal new results are those for dual-lift and 
multilift systems. Equations for three dual-lift config- 
urations are given and it is shown that these can be 
integrated in a single simulation. The multilift sys- 
tem extended to any number of helicopters was al- 
sotreated. These previously difficult problems are seen 
to be tractable for derivation, analysis, and program- 
ming by hand with the methods of this report. 


Equations for degenerate-body approximations 
(point masses, rigid rods) are given to accommodate 
various simplifications useful in control analysis, and 
results for dual- lift and multilift systems are given. A 
reduced-order load-suspension system is obtained by 
assuming that the helicopters control the system solely 
through their applied forces with the helicopter atti- 
tudes in steady state relative to the reduced system 
motion. The reduced system equations are those ob- 
tained assuming point-mass helicopters. The reduced 
system can be studied and a control law formulated 
for the applied forces independently of any helicopter 
details and of the problem of implementing the applied 
force controller in the helicopters. This simplification 
is expected to be realistic for helicopters with single- 
point attachments, including the dual-lift and multilift 
systems. 

Linearized equations of motion for general slung- 
load systems with inelastic suspensions are derived 
from the nonlinear equations of motion, and results 
for the dual-lift system with spreader bar are given for 
general static-equilibrium flight conditions. 
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APPENDIX A 

SUMMARY OF USEFUL KINEMATIC RELATIONS 


INTRODUCTION 

This appendix collects the general kinematic formu- 
las used in applying the methods of this report to de- 
rive the simulation equations given in appendixes B-E 
and the linearized equations for dual-lift systems. 
These formulas are mostly counterparts of relations 
from the classic physical vector theory of rigid-body 
dynamics given here in a form appropriate for digital 
computations by introducing coordinate-frame trans- 
formations. They facilitate derivation of the simula- 
tion equations in the preferred form of this work, in 
which all terms appear as operations on the underly- 
ing physical vectors. 


dV i 

nt 01 




time-derivative of physical 


vector V relative to T\ and given by the 
coordinates in T\ 


The notational specializations stated in the list of sym- 
bols are continued in this appendix (e.g., Va, u b are 
reserved for the inertial velocity of Ra, j F bl respec- 
tively), but there is no conflict with the notation shown 
above, in which superscripts are used. Definitions and 
specializations reflect aeronautical usage and are gen- 
erally consistent with those in reference 41. 


Skew-Symmetric Matrix 


The notation for this work is stated in the list of 
symbols. Occasionally, the general notational rules 
for vector kinematics advocated by Kane (ref. 40), are 
used; these are as follows: 

Velocity of point Ra relative to reference 
frame T\ 

Angular velocity of reference frame T\ or 
of rigid body B\ relative to frame Ti 

Time-derivative of physical vector V rela- 
tive to reference frame, T\ 

To work with the scalar equations used in digital com- 
putations a subscript can be added to indicate the ref- 
erence frame in which the coordinates of the vector are 
given: 

1 Vj 2 Velocity of point Ra relative to frame T\ 

given by its coordinates in frame JF 2 

2 o >3 Angular velocity of frame T\ relative to 

frame Tz given by its coordinates in frame 

Time-derivatives of these objects are necessarily with 
respect to the reference frame in which they arc given: 


lya 


2 u , 1 


dt 


First, the general skew-symmetric matrix S(x, y, z) 
is defined from the scalar triplet {x,y,z) as shown in 
table 6. This allows scalar representation of the vector 
cross products, which occur frequently in this work, as 
shown in equation (51) for vectors VI, V2 referenced 
to frame T a > The algebra of skew-symmetric matrices 
is consistent with corresponding relations in vector al- 
gebra such as product reversals (eq. (51)) and triple 
products. Geometrically, S(Vl a ) maps any vector V 
given in T a to the vector VI x V referenced to T a 
and perpendicular to the plane of (VI, V). The refer- 
ence frames of vectors occurring in expressions based 
on equation (51) can be selected in any convenient 
and consistent way using transformations (eq. (52)). 
When viewed as an isolated matrix, the columns of 
5(14) physically represent cross products of V with 
the axes of T a and referred to T a (eq. (53)). Identities 
in the matrix 5(14) represent vector operations on an 
arbitrary vector, such as product reversal and the cross 
product with a sum of vectors (eqs. (54) and (55)). 

Cross products representing Coriolis velocities and 
accelerations and centrifugal accelerations (eqs. (56)- 
(58)) are the basis from three-dimensional kinematics 
of virtually all terms in An, A~ l v, Au y in the appli- 
cations of this report. Cross products also represent 
the applied moments due to cables in the term /c. In 
general, the force F applied to body B at point Ra 
imposes a moment on B about the point Rc, which is 
given by = Rea x F. In the present work, this 
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can be specialized to cables in tension and to moments 
about the rigid body c.g.’s (eq. (59) and accompanying 
sketch). If cable Cj applies tension TCj in the direc- 
tion kcj at point Rj on body Bi , then its moment 
about the c.g. of Bi is given by 

Mlj = Ri*j x kcj TCj = £ij TCj 

where the symbol £ij is reserved for the moment of Cj 
on Bi per unit tension. 


Transformation rates arise in the inertial coupling 
term A or A u in the equations of motion, and useful 
expressions are noted in equations (69) and (70). In 
general, if T ayb is a transformation from T b to T a , then 

T a , b = T 0l 6 (79) 

This relation follows as a counterpart of the Coriolis 
equation relating the time-derivative of a vector, V, 
relative to two different frames, J- a , Tb\ that is, using 
Kane’s notation, 


Coordinate Transformations, Angular Veloci- 
ties, and Transformation Rates 


a dV 

dt 


b dV 

dt 


+ a u> b x V 


Formulas defining the transformation of a vector 
from its coordinates in inertial space, 1 to its co- 
ordinates in body axes, are given in table 7 along 
with relations governing its occurrence in kinematic 
equations with time-derivatives and its gradient. 


but the time- derivative of 

V a = T at b V b 


is 

V a = T a ,b Vb + T Qi b Vb 


The primitive transformations for rotations about 
a single axis i,j,k of a right-handed orthogonal ref- 
erence frame, T are denoted /^(tf), ^ 3(^)1 re ~ 

spectively (eqs. (60)- (62)). Then the usual Euler-angle 
transformations of aeronautics is given by the yaw, 
pitch, and roll sequence of rotations T^jv as illustrated 
and defined in table 7 (cq. (63)). It is often useful to 
note that the rows of T b ,N are the Ts -components 
of the axis vectors of T b and its columns are the T b - 
components of the axis vectors of Tn (eq. (64)). The 
notation (0,0,^) is reserved in this work for inertial 
angles, and a subscript is attached to indicate the body 
axes. Euler-angle transformations from other axes T a 
are sometimes needed and an ad hoc notation is de- 
fined in text (e.g., A <f>b, A 0^ A V^)- 

In general, if frame T b is obtained from Tn by 
any arbitrary sequence of independent rotation angles 
/3 U /? 2 , . - . , 0n about the axes ul, u2, . . . , un, respec- 
tively, then its inertial angular velocity is given by the 
superposition 


whence T 0j t> Vb is identified as the scalar counterpart of 
a u; b x V for arbitrary V b \ this establishes equation (79). 
In all cases, A u is composed of terms that are Cori- 
olis and centrifugal accelerations that arise from the 
Coriolis equation and utilize equation (79). 


Derivatives of transformations with respect to Euler 
angles occur in the linearized EOMs, especially in the 
position-perturbation term. These derivatives can be 
obtained from equation (79) by identifying the terms 
in the expansion 


Tb^N 


QTb, N ; QT byN 

94>b 06 dOb 


0 b 4* 


&Tb,N 

d'lpb 


Tpb 


with corresponding terms in the expansion of 
Tb y N S(vbw), using u?b = <£& ib + 9 b jb' + <j) b kN and 
the distributive property (eq. (55)). The results are 
listed in table 7 (eqs. (71)-(73)) and here as 


dTb t N 

d<j) b 


= -T ht N S(ib N ) = - S{ib b ) T btN 


N u; b = $1 ul + • ■ . , +/?6 un (78) 

Equation (66) specializes this rule to the Euler-angle 
sequence defined above, and the familiar linear relation 
from aeronautics for u;bb{c*b) and its inverse arc noted 
in the table, along with expressions for Euler-angle 
rates in terms of dot products with u>b. The notation 
Wb is introduced to indicate the row list of the axes 
of roll, pitch, and yaw rotations for body axes T b , and 
a subscript indicates the coordinate frame in which 
these arc given. 


dT bi N 

d0 b 


= -T b ,N S(jb' N ) = -5(j6j,) T btN 


(80) 


-zrr- = -Tb t N S(kN N ) = -S(kNb) T b , N 

oi> b 

These derivatives are cross products with the axes of 
roll, pitch, and yaw rotations, which are illustrated in 
table 7. They are also maps from Tn to T b , and the 
rotation axes can be given in any convenient reference 
frame by using transformations as required for consis- 
tency with the input and output reference frames. 
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Derivatives of transformations arise commonly in the 
linearized EOMs from gradients of vectors, 


V b = V N , V N = T N ,b Vb 

relative to Euler angles, where the vectors on the right- 
hand side are independent of #6, It follows 
from equation (80) that 


OTb,N VN 
d<j> b 


= -(ib x v) b = (v x ib) 6 


dTb,N vn 
d0 b 


= ~(jb' x v) b = (v x jb')6 


(81) 


dT b ,N V N 
dipb 


— -(kN x v)b = (v x kN)b 


The vector-mechanical meaning and output reference 
frame is apparent in equation (81). Reference frames 
for the vectors and the order of the cross product can 
be selected in any convenient way consistent with these 
characteristics. By using equation (81), the gradient 
with respect to the Euler-angle triplet is 


Vab(Tb,N v n) = T by N S(v N ) [ib N ,jb' N ,kN N ] 


= T bt N S(v N ) Wb N (82) 


and also 


VL(Tjv f 6 t>6) = -Tn, b S{v b ) Wb„ = -S(v N ) Wb N 


Again, the notation Wbn,Wb b indicates the matrix of 
Euler-angle rotation axes given by inertial or body-axis 
components, respectively. These results arc included 
in table 7 and are valid for any arbitrary vector, v. 

Cables Axes 

The general treatment of cable angles and cable axes 
in this report is summarized in table 8. The cable di- 
rection kc is located by inertial roll and pitch angles 
((f) c 9 C ) taken in the usual Euler sequence, as illus- 
trated in the table and in equation (83). 


Cable axes T c = {ic,jc,kc} are constructed from 
these two angles with kc along the cable and ic in the 
inertial vertical plane of (iN,kN). The transforma- 
tion T Cy N (eqs. (84) and (85)) is a specialization of the 
usual Euler-angle transformation with i}) c = 0, and its 
rows and columns are the axes of T c expressed by its 
coordinates in Tn and, conversely, as before. 

Useful formulas are given for the inertial angular ve- 
locity of T c in terms of cable-angle rates (eqs. (86)- 
(88)), and for the cable velocity expressed in terms of 
cable length and angle rates (eqs. (89)- (92)). The T c 
coordinates of the cable velocity (eq. (92)) separate 
the cable-stretching motion (l c ) from the orthogonal 
motion caused by cable rotation; this fact is frequently 
used in selecting the generalized coordinates u in the 
applications. 

Cable angles relative to some noninertial axes T\ 
(e.g., helicopter body axes or level-heading axes) may 
be more useful coordinates in some problems. For 
these cases, Euler pitch and roll angles relative to 
T\ are used with special notation (e.g., A (f) c , A0 C ) 
to distinguish them from inertial angles. The previ- 
ous coordinate transformations and illustration apply 
by analogy with a change of notation (eq. (93)). The 
longitudinal cable axis ic is now in the T\ vertical 
plane of (il,kl). The velocities of Rc, T c relative 
to T\ (eqs. (95) and (96)) are given by analogy to 
equations (88) and (92). The relative cable velocity 
in cable axes l V* separates the cable-stretching mo- 
tion and motion due to cable rotations relative to T\ 
into orthogonal components. The inertial transforma- 
tion and velocities (eqs. (97)-(99)) are obtained from 
elementary rules. Cable-stretching motion is again iso- 
lated in Vc and its rotational motion is represented a s 
a superposition of orthogonal rotation relative to T\ 
plus the effect of .TVs inertial angular velocity. 

Gradients of the cable-axis transformations T Ci n oc- 
cur in the linearized EOMs. The required general for- 
mulas are specializations of the previous results ob- 
tained by dropping the yaw derivative and the axis of 
yaw rotation kN from table 7 (eqs. (71)-(75)). 
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Table 6. Skew-symmetric matrices and cross products 


i 

1 

o 

V \ 

II 

Co 

2 0 

- -x 

1 

V -y x 

0 / 


(VI x V2)„ = S(V1 0 ) V2 a = — (V2 x VI). = -S(V 2 a ) Vl a 
(VI x V2) a = S(V1 0 ) T a , b V2 b = T a<b S(Vl b ) V2 b 
S(V a ) = [(V x ia) a , (V x ja) ai (V x ka) a ] 

S T (V a ) = -S(V a ) 

5(V1„ + V2 0 ) = 5(VT a ) + S(V2 a ) 

Coriolis and centrifugal terms: 

(a; x R) a = •S'(u' a ) R a = —S(R a ) u> a 
(a; x V) a - S(ui a ) V a = -S(V a ) w a 

(w x a; x R) a = 5 2 (w a ) R a = -5(w a ) S(R a ) u> a = u> • Rw a — w 2 i? a 

Moment of cable Cj on body Bi about Ri* : 

Ti = {ii, ji, ki} = body axes for Bi 

Miji = (Ri‘j x TCj kcj)i = S(Ri*ji) kcj, TCj = £ij t TCj 


(50) 

(51) 

(52) 

(53) 

(54) 

(55) 


(56) 

(57) 

(58) 


(59) 
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Table 7. Coordinate transformations, angular velocities, and transformation rates 


Single-axis transformations: 


no ox 


Ex(a) = 

0 COS <7 

sin 

<7 

( 60 ) 


\0 — sin a 

COS 

a) 



/cos a 

0 

— sin 



E 2 (<t) = 

0 

1 

0 


( 61 ) 


V sin a 

0 

COS (7 / 



/ COS (7 

sin <7 



E 3 (ct) = 

— sin a 

cos a 

0 

( 62 ) 


K 0 


0 

J 



Euler-angle transformation, 

T n = {iNJN.kN} = inertial axes 

Ti — {ib,jb,kb} = body axes 

ab — Ob, ip b) T = Euler-angle triplet 


kN 



(a) = heading rotation 

(about kN) 


ib 



(b) = pitch rotation 

(about jb") 


ib 



(c) ^ = roll rotation 
(about ib) 


Euler angles. 
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Table 7. Continued. 


T biN (ab) = Ei(4 > b ) E 2 {0 b ) Ei{ip b ) 




cos ifib cos 6b 

sin 4>b cos ifib sin 06 — cos (fib sin ifib 
cos (fib cos ifib sin 6b 4- sin (fib sin ifib 


sin ifib cos #6 

sin (fib sin ifib sin #5 -f cos (fib cos 06 
cos 06 sin 06 sin 9b — sin 06 cos ifib 


— sin 06 \ 

sin (fib cos 
cos (fib cos / 


( ) 

Tb,N “ 

Vfcbjy / 
^iV,6 = 7^ 


= (W 6) jAT 6 , fciV 6 ) 


Angular velocity of Tb relative to 
N w b = i> b kN + e b jb' + <j> b ib 



f 1 

0 

— sin 6 b \ 

/M 


/M 

uibb = Wbb cxb = 

0 

COS (fib 

sin (fib cos Ob 

@b 

= (ib b , jb' b , kN b ) 

o b 


Vo 

— sin (fib 

COS (fib cos 06 / 

\ipb/ 


\i>b) 



z 1 

sin <fib tan 6b 

cos 06 tan 6*6 \ 


/ ib 7 • bo b/ cos 06 \ 

ab = W^ 1 w& 6 = 

0 

COS (fib 

— sin (fib 

= 

jb' • a;b 


Vo 

sin (fib / cos 6b 

COS 06 / COS / 


\ kb 7 • (Jbf cos 06 / 


where {ib'Jb',kb'} are unit vectors related to J 7 ^ as illustrated above. 


Transformation rates: 


T a ,b = T a>b S(“u;£) = S( a u> b a ) T atb 


TN,b = TjV.t S(ub b ) = S(ujb^) Twfi 
T b ,N = -S(u>b b ) T bt pj = —T bt N S(ojbN) 


(63) 

(64) 

(65) 

( 66 ) 

(67) 

( 68 ) 


(69) 

(70) 
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Table 7. Concluded. 


Transformation gradients: 

= (^) T = -n,N s(ib N ) = -sw 

W = (^if) T = S(jb' N ) = -5(j6') r 6)JV 

= (%^) T = —Tb,N S(kN N ) = —S(kNb) T b>N 

V£ b Tb'N(ub) vn — Tb,N S{vn) W frjv — S(vb) W 6t 
T N , fc (a6) v b = -T h r, b S(i>b) Wb b = -S(vjv) t'F6 w 


Wb = [ib, jb',kN] = roll, pitch, yaw axes of rotation 


Wb N = 


Wb b = 


/COS ^6 COS Tpb 

cosOb sin ^6 
\ — sin #£, 

1 0 

0 cos 4>b 
VO — sin <pb 


— sin ipb 0 \ 

cos xpb 0 

0 1 / 

— sin \ 

sin 4>b cos 6/j 
cos (pb cos 6 b ) 


(71) 

(72) 

(73) 

(74) 

(75) 


(76) 


(77) 
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Table 8. Inertial cable angles and cable axes 


Cable direction: 

N = {iN, jN, kN} = inertial axes 
{@ c , 4>c) = cable pitch and roll angles 
kc = cos <p c sin 9 C iN - sin 0 C jN + cos 4> c cos 9 C kN 


ic 



(a) Pitch rotation (about jN) 


(83) 

ic 

jN 
jc 

kc' 

(b) Roll rotation (about ic) 



Cable angles. 


Cable axes, F c \ 


F c = {ic, jc, kc} 


T Ct N = E 1 (d> c ) E 2 (0 c ) = 


( cos 9 C 0 — sin 0 C \ 

sin 4> c sin 0 C cos (fi c sin 4> c cos Q c 

V cos 4> c sin 9 C — sin 0 C cos (f> c cos 9 C / 


L c,N 


( ic N \ 
JCN 

\kc T N J 


(iN c jN c kN c ) 


(84) 


(85) 
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Table 8. Concluded. 


Inertial angular velocity of T c : 

u;c = 9 C jN + 0c ic 

ujc = 0 C ic -f 0 C (cos 0c jc — sin 0 C kc) 

o;c = 0 C (cos 0 C iN — sin 6 C kN) -f 0 C jN 

Inertial velocity of cable line segment, Rc: 


Rc = i c kc 

Vc = i c kc + t c kc 

kc = ojc x kc = 0 C cos 0c ic — 0 C jc 

Vc = i c e c COS 0c ic — 0c l c jc + i c kc 


Cable angles relative to non-inertial axes, T\: 


= {iljl.kl} 

(A9 C , A<p c ) = cable pitch and roll angle relative to T\ 


T c>1 = Ei(A<pc) E 2 { A9 c ) = 


{ cos A 9 C 
sin A9 C sin A<f> c 
V sin A6 C cos A<p c 


0 — sin A6 C \ 

cos A<f> c cos A 9 C sin A <p c 

— sin A<p c cos A9 C cos A 4> c J 


ic = cos A 6 C il — sin A 9 C kl 

•*-oi c = A 9 C jl + Ad> c ic = A<fi c ic + A 9 C (cos A <j> c jc — sin A4> c kc) 

= A9 C jl + A<p c (cos A 9 C il — sin A 6 C kl) 
V c = £ c kc + d c ui c x kc = (£ c A 9 C cos A <p c ic — i c A0 C jc -f l c kc) 


T c ,n = T c> i Ti'W 

UJC = + wl 

Vc = 1 V c + l c Ul x kc 


( 86 ) 

(87) 

( 88 ) 


(89) 

(90) 

(91) 

(92) 


(93) 

(94) 

(95) 

(96) 

(97) 

(98) 

(99) 
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APPENDIX B 

SIMULATION EQUATIONS FOR SINGLE-POINT SUSPENSIONS 


INTRODUCTION 

Figure 4 shows three suspensions of interest with a 
single attachment point on the helicopter. If the cables 
are modeled as inelastic, then the multicable suspen- 
sion with three or more cables (fig. 4(a)) imposes three 
constraints on the load motion whereas the other two 
suspensions impose only one. If cable elasticity is mod- 
eled, then the suspension in figure 4(b) is considered 
to have only one elastic cable; the sling legs that con- 
nect the load to the ring are assumed to be sufficiently 
short and stiff that their elastic stretching is negligible. 
Figure 4(c) represents suspensions with long sling legs 
whose elastic stretching must be considered. 

Previous simulations of such systems are described 
in references 13, 14, and 16. Reference 14 documents 
a Langley Research Center simulation of the CH-54 
helicopter with a MILVAN cargo container suspended 
from a single elastic cable as in figure 4(b). Refer- 
ence 13 uses a general formulation for elastic suspen- 
sions with multiple attachment points in which cveiy 
cable connects the two bodies; this is readily special- 
ized to systems a and b in figure 4, and an approximate 
adaptation to system c, in which no cable connects two 
bodies, is given. Reference 15 considers a system like 
a with elastic suspension and one like system b with 
elastic or inelastic suspension. These formulations all 
begin with the rigid-body accelerations, either inertial 
or relative to body axes, and relative load velocity is 
calculated in some cases. 

This appendix contains simulation equations for sys- 
tems a and b obtained by the methods presented in this 
report. Generalized velocity coordinates are selected 
specific to each case in order to separate the system 
motion caused by cable stretching from motion with 
invariant cable lengths. These are rigid-body veloci- 
ties and cable velocity or relative motion coordinates in 
an appropriate coordinate frame. The results account 
for both elastic and inelastic suspensions. Interaction 
forces are explicitly calculated in both cases. Appro- 
priate partitioning coordinates for system c were not 


obtained, but coordinates that simplify the determi- 
nation of the suspension forces for elastic cables are 
indicated without elaboration to EOMs. 

Nomenclature and enumeration of the attachment 
points, cables, rigid bodies, and system parameters 
used in this appendix are defined in figure 4. The 
cable-length parameters {£oj} refer to the fixed ca- 
ble lengths in the case of inelastic cables and to the 
unloaded cable lengths in the case of elastic cables. 
Controllable parameters are not considered (p = 0). 


Multicable Suspension 

System and constraints- Referring to figure (4a), 
77 i cables, (71,(72, . . . , Cm attach the load, 52 to the 
point Ra on the helicopter, 5 1. Usually, three or four 
cables are used. In any case, it is assumed that the 
cable geometry is such that if the cables are inelastic 
then the suspension fixes the distances from Ra to 
three noncolinear points on 52. This suffices to impose 
three holonomic constraints on the configuration 

c = 3, d = 9 (100) 

by fixing the line segment between the load c.g. and 
the helicopter attachment point in load-body axes; 
that is, for inelastic cables, 

i?a2* 2 = 0 

where 

Ra2* 2 {r) = T 2jN (R^n - fil V ~ t n . i 

For some arrangements of three cables, collapse of a 
cable is possible, thus leaving the load free to rotate 
about a line joining the remaining two attachment 
points, in which event the number of constraints is 
reduced to two. However, this is prevented by adding 
a cable (m = 4), and then c = 3 even if one cable 
collapses. With four or more inelastic cables the con- 
straints are imposed redundantly, and the individual 
cable forces cannot be calculated; however, their resul- 
tant can always be calculated, and this calculation is 
independent of the number of cables used to maintain 
c= 3. 
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Parameters: 

ml, Jl, Ri* ai 
m2, J2, R2*l2, R2*m2 
j "I, m for m > 2 

(a) Multi-cable 
suspension 



ml, Jl, Ri* ai 
m2, J2, R2*1 2 
Ho 

(b) Single-cable 
suspension 

Figure 4. Single point suspensions. 


ml, Jl, Rl* ai 

m2, J2, R2*2 2 , ..., R2*m 2 

&>j> j=1, m for m > 2 

(c) Single-cable sling 
suspension 



Generalized velocity coordinates and config- 
uration kinematics- The first task is to find gener- 
alized velocity coordinates for the elastic system 



such that, if the cables are inelastic, then A = 0 and 
ul are the generalized coordinates. In the present 
case, (V 1 * n ,u) 1 uU j 22 ) comprise nine generalized ve- 
locity coordinates which suffice to define the configu- 
ration velocity of the inelastic cable system and can 
be taken as ul. The remaining element of the con- 

guration velocity, V2* N , is given from the inertial 
derivative of 

R2 N = Rl* N + t na RVai + T n<2 Ra2* 2 (102) 


as 


' - Tn,i S(Rl*ai) wli 

~T n ,2 S(Ra2* 2 ) uj2 2 
(inelastic cables) 


VZ*N={ 


*n~Tn,i S{RVa x ) oilj 
+Tn , 2 [Ra2\ - S(Ra2* 2 ) w2 2 ] 

• (elastic cables) 

T , , „ (103) 

inus, the configuration velocity caused by cable 

stretch can be given by augmenting the coordinates 
u for the inelastic cable system given above with the 
velocity of Ra2* relative to load-body axes, Ty 
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/™V\ 


u — 


U) 2 2 


\Ra2\J 


(104) 


The inverse relation for Ra2 2 (v) is obtained from equa- 
tion (103) as 

Ra2* 2 = T 2iN (- Vl* N + F2* n + T^i S(fll*ai) wh) 
-\-S(Ra2* 2 ) v2 2 (105) 


The complete kinematic relations v(u ), u(v) are as- 
sembled in the equation summary, table 9, using equa- 
tions (103) and (105). The submatrices All , A of A~ l 
required below are identified in the table. Only coor- 
dinate transformations and skew-symmetric matrices 
representing Coriolis terms occur in A, A -1 . Other- 
wise, these matrices each contain nine rows from the 
unit matrix corresponding to the nine coordinates of 
u, which are also configuration velocity coordinates. 

External forces and inertia coupling terms- 
The external forces and inertia coupling terms fo 
are assembled in part b of table 9. The notation 
F01,...,M02 for the vector elements of fo is in- 
troduced for brevity in later equations; fo contains 
the sum of the applied forces and moments due to 
weight (ml g, m2 g), aerodynamics and rotor out- 
put (FA1, . . . , MA2), and the inertia coupling terms, 
X + D A u. The time-derivative of A requires only 
time-derivatives of transformations, for which a gen- 
eral formula is given in appendix A. The three terms 
in F02^ due to Au are recognized as centrifugal and 
Coriolis acceleration terms of the form wxwxR and 
u; x V. 


Suspension forces- The forces applied to the con- 
figuration of rigid bodies by the suspension are denoted 
for this system by 


(FC l N \ 


fc = 


FC2 n 

MCI i 


\mC2 2 J 


(106) 


where FC1, . . . , MC2 are the resultants of cable forces 
acting on each rigid body and their moments about its 
c.g. 


For inelastic cables, fc is given by equation (17) as 


/ ~Tn, 2 \ 


/c = As = 


Tv, 2 


— A22 T/v, 2 



V — S(i?a2* 2 ) / 


(107) 


where A T are the last three rows of A -1 and s are 
suspension force parameters. Since FC In = —Tv, 2 s , 
5 is identified as the T 2 - components of the resultant 
suspension force applied to the helicopter at Fta; that 
is, s = — FC12- In that case, s can be conveniently 
replaced by -T 2 ,v FCIn in equation (107), and T 2 ,v 
can be combined with A to obtain the simple, alterna- 
tive expression, 


l 1 \ 


fc = A' FCIn 


-I 

a 22 


FC 1 N 


\AlJ 


(108) 


Finally, for inelastic cables, FCIn is given in terms of 
A' and fo from equation (19). The results are listed 
in the equation summary, part c of table 9. 


For elastic cables, fc can be given by equation (11), 
which applies generally whenever each cable connects 
two rigid bodies. However, it is simpler to use the 
form in equation (108), which is valid in this exam- 
ple whether or not the cables are inelastic. It only 
remains to calculate the resultant suspension force, 
FCIn, from the cable tensions: 

FCl N = kcj N TCj 

= £7=1 ™ a *{°’ K i (1 - T?) + c, %}RajN 

(109) 


Cable damping, ? is usually absent from slung-load 
simulation models. The cable line segments required 
in equation (109) can be calculated from 

RclJn = Tv, 2 (RaT 2 + R2* j 2 ) j = 1 , ..., m 

and each cable length rate can be obtained as the com- 
ponent of cable velocity along the cable, and can be 
shown to be 

tj = RajJ Ra2* 2 jtj j = 1, ...,m 

where {R2* j 2 , j = 1, . . . , m) are given system parame- 
ters and where Ra 2* 2 , Ra2\, 7\-,2 are obtained directly 
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from the generalized position and velocity coordinates. 
These results are included in the simulation equation 
summary (table 9). 

Simulation equations- Finally, the simulation 
equations are listed in part d of table 9. The total 
specific force, sf, due to fo + /c, is assembled and its 
vector elements are denoted SF1,,..,SM2 for conve- 
nience. Last, the elements of u are obtained by ex- 
panding A~ l sf. If the cables are inelastic it is unnec- 
essary to calculate Ra2* 2 or the term due to Ra2* 2 in 
A it, but they can be evaluated to monitor computa- 
tional accuracy. 


Single-Cable Suspension 

System and constraints- The suspension shown 
in figure 4(b) has a single main cable, C, attached to 
the helicopter at Ra. Additional cables or hardware 
attach the load to Rl, but these are assumed inelastic 
here, being either short cables with negligible stretch, 
or inelastic links of any length. In these cases, they 
can be regarded as part of the load rigid body, B2, 
since the point Rl is fixed relative to the load. 

If Cl is inelastic then it imposes a single holonomic 
constraint on the configuration motion by fixing the 
distance between Ra on B 1 and Rl on Z?2, in which 
case 

c= 1, d= 11 (110) 

and the constraint equation can be given as 

£(r) = |Ral*r| 

= \R2*N + T N ,2R2*l 2 -Rl\-T Ntl JH*a! | 
= to 

Generalized coordinates- The first task is to 
select appropriate generalized coordinates (i/l, A) 
with the properties previously stated with equa- 
tion (101). Assume that the nine rigid-body coordi- 
nates, (Fl* N ,u;li,u> 22 ) can again be included in id. 
Then the remaining rigid-body velocity, 1^2*^, can be 
given from the derivative of 

7^2*^ — -b Tjv,i Rl*a\ 4- Raljy — T)v,2 R%* I 2 

as (111) 

V2*n = ^1%/ ~~ Tn,i S(Rl*ai) u^l 1 + Val^ 

+Tv,2 S(R2*1 2 ) u>2 2 


It remains to express Val N in terms of coordinates 
that separate the cable rotational motion from the ca- 
ble stretching motion, £. The rotational coordinates 
can be defined in several ways that represent rotations 
relative to either inertial space or to helicopter 
body axes T\, or to level heading axes based on the 
helicopter heading Th - The kinematics are given for 
all three cases, and the simulation equations are given 
for coordinates defining cable rotations relative to 

First, inertial cable pitch and roll angles 9 c ,<p c 
are defined in table 10, and then cable axes Tc - 
{ic,jc,kc} can be constructed with kc along the ca- 
ble and ic in the inertial vertical plane of (iN,kN). 
The angular velocity of T c relative to J~n is, from ap- 
pendix A; 

ojc = 4> c ic -f 0 C (cos 4> c jc — sin 0 C kc) 
and the velocity Val is given from 

Val = £ kc + £ ujc x kc 
or 

Val^ = T N , e Val c — T^ iC (£ & c cos0 c , — £ <j> c , £) T 

( 112 ) 

The coordinates Va l c are seen to separate the cable 
rotational and stretching motion and, together with 
(Vl* N ,wluv2 2 ), they are suitable generalized coordi- 
nates for the elastic cable system. The cable rotation 
is represented by linear velocity components perpen- 
dicular to the cable direction. This results in a sim- 
pler kinematic relation, v(u] } than is obtained by us- 
ing the cable angular rates 9 C ,4> C as coordinates. The 
kinematic relations u(w), it(n) are readily given from 
equations (111) and (112), and the results are listed in 
part a of table 10. 

Most loads suspended by a single cable will stabilize 
at a steady-state position relative to the helicopter on 
each steady segment of a reference flight path. For 
such loads, in a steady turn, the cable will trail the ver- 
tical at a steady angle in the vertical plane of (il, kN) 
because of load drag, and the load will swing out to a 
steady cable angle that is about equal to the helicopter 
roll angle, <f> 1 , from this vertical plane as a result of cen- 
trifugal force. These steady angles are represented by 
sinusoidal variations of the inertial cable angles 0 C , 0 C 
with helicopter heading, ipi. For example, 9 C alter- 
nately represents the trailing angle or the swing angle 
at different headings. To avoid this complexity in rep- 
resenting ordinary steady motion, other choices of the 
generalized coordinates are considered next. 
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Second, cable pitch and roll angles relative to heli- 
copter body axes A 9 C , A 0 C , are defined in part b of ta- 
ble 10. The following derivation uses general formulas 
for cable kinematics relative to noninertial axes given 
in appendix A, table 8. Cable axes, T c — {ic,jc,kc} 
are again constructed with kc along the cable, but ic 
is now a direction in the helicopter vertical plane of 
(il,kl). In this case, a steady load-trail angle is rep- 
resented by a steady relative pitch A0 C , and a steady 
load-swing angle is represented by a relative roll an- 
gle that is approximately zero. For these coordinates, 
the angular velocity of T c relative to T\ is (table 8, 
eq. (95)) 

Awe = A0 C ic + A# c (cos A <f> c jc — sin A<f> c kc) 
and then the inertial velocity Val is (table 8, eq. (96)) 
Val = e kc + e(u>l + Awe) * kc 
or 

Fnljv = — ^Nyi S^Rali) uj li T AFnl c 

where 

AVal c = {£ A0 C cos A^ c , -l A<£ C) i) T (113) 

Here, AVal is the velocity of Ral relative to 
T\, and its components in T c are seen to separate 
the cable rotational and stretching motions; then 
(^i v uli, u2 2 , Altaic) are suitable generalized ve- 
locity coordinates for the elastic system. The kine- 
matic relations t?(u),u(v) are given by using equa- 
tions (111) and (113), and the results are listed in 
part b of table 10. 

Third, level heading axes, Th = {ih,jh, kN}, are a 
local vertical frame defined from the helicopter heading 
ipi, as shown in part c of table 10. Cable pitch and roll 
angles relative to Th, are denoted 9 c h,<j> c h, and then 
cable axes T c — {ic,jc,kc} are constructed (part c 
of table 10) with ic in the inertial vertical plane of 
(ih,kN). In this case, 9 c h,4>ch have steady values in a 
steady turn that are approximately the load- trail and 
load-swing angles, respectively. The angular velocities 
of Th relative to Tn, and of T c relative to Th are 

uj h = ip i kN 

uch = 4>ch ic -b 6 c k{ cos 9 ch jc - sin <p c h kc) 

where ip\ is given in terms of o>li by the usual kine- 
matic relation for Euler-angle rates (eq. (68) in ap- 
pendix A): 

Ipi = — l — (0, sin^i, cos^i) wli = — b- (fcli) T wli 

COStq COS^I 


where kl' can be shown to be a direction in the vertical 
plane containing (il,ih,kN) at an angle 9\ from kN. 
Then the inertial cable velocity is 

Val = i kc -b £ (u;h + ejeh) x kc 

or 

Faljv = Tjvc ( — “T (kN x kc) c (kl\) T ujli+Valh c ) 

’ cos 9 1 

(114) 

where 

Vd\h c — (-^ Och COS (pchi ^ 4*ch, 0 

Here, Valh is the velocity of Ral relative to Th, 
and its components in T c are seen to separate the ca- 
ble rotational and stretching motion. The coordinates 
(Vl* N ,a;li,u;22, Valh c ) are suitable generalized coor- 
dinates for the system with elastic cables. The kine- 
matic relations v(u),u(v) are readily given by using 
equations (111) and (114), and the results are listed in 
part c of table 10. 

The remaining equations for simulating the system 
with these coordinates are given in parts d-f of ta- 
ble 10. The suspension force is obtained from equa- 
tion (11) using the last row of A -1 which corresponds 
to £(v). 


Remarks 

1. In the event of extraordinary motion in which all 
the cables collapse (fc — 0), the equations for elas- 
tic cables still correctly represent the motion of the 
two independent rigid bodies. However, the equations 
for inelastic cables assumed c constraints. If the num- 
ber of constraints is reduced as a result of collapsed 
cables, the d equations no longer suffice, but the simu- 
lation can accommodate this regime by carrying along 
the complete set of 12 equations. However, additional 
equations are needed to detect the onset of such col- 
lapse, but these are not given here. Since such extreme 
divergent load motion is unacceptable it may be un- 
necessary to simulate it. 

2. The suspension of figure 4(c) has not yet been 
specifically discussed. This suspension consists of a 
three- or four-cable sling which attaches the load to 
the single cable Cl, and the elasticity of the sling is 
to be simulated. If the sling cables are inelastic, then 
the treatment and equations of table 10 apply, since 
the load and sling still form a rigid body. If the sling 
cables are elastic, then all cables must be considered 
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in order to determine fc. In this latter case, each ca- 
ble is connected to only one of the rigid bodies. The 
principal analytical problem is to locate the intercon- 
nection point, Rl, from the force-balance equation at 
R1 and from the sling geometry in order to obtain 
cable lengths and tensions: 

FC1 2 = TCI kcl 2 = K x (1 - -~)Ral 2 

= ~ T Kj (1 - —■)} Rjl2 

Ra\ 2 = Ra2 m 2 + R2*l 2 = nkcl 2 (115) 

Rjl 2 = R2*h - R2* j 2 = £j kcj 2 j = 2 , . . . , m 

Cable Cl can be taken as elastic or inelastic in the 
force-balance equation, and cable damping is omitted. 
If Ra2* 2 is known as a result of using the coordinates 
v or including Rci2* 2 In the generalized velocity co- 
ordinates, then equations (115) is a nonlinear vector 
equation in R2*l 2 and the (dependent) cable force, 
FC1 2 . Reference 13 provides an approximate solution 
for R2*l 2 . 

Alternatively, if the generalized coordinates u are se- 
lected as (VT*^, odi, uj2\ , R2*\ 2 )^ then R2*\ 2 and all 
sling-leg lengths and tensions are known, and FC\ 2y 
Ral 2 are readily given from equations (115). This use 
of specially selected generalized coordinates for the 
elastic-sling case circumvents the problem of solving 
equation (115) for R2*l 2 . Moderate complexity reap- 
pears in the velocity relations V2* N (u) and R2*l 2 (v) ) 
which now depend partly on the sling geometry and 
elasticity parameters: 

V2* n = W N - T Ni i S(RV ai ) u>h 

^T/v,2 S(Ra2* 2 )cj2 2 — T/v,2 [f T M] R2*l 2 

or, solving for R2*l 2 , 

R2*l 2 = -[/ + A/]" 1 T 2 , n [V2* n - VV N 

+Tjv, i S(Rl*ci\) lj\\ + Tjv , 2 S(Ra2* 2 ) u;2 2 ] 

(116) 

where 

Ra2* 2 = Ra\ 2 - R2*l 2 

and M gives Ra\ 2 from R2*\ 2 . Noting in equa- 
tion (115) that Ral 2 can be given from FC \ 2 which, 
in turn, can be given from R2*l 2 , then 


Ral 2 — Ral 2 ] FC1 2 

= [Vr ci3 flal 2 ][V^ 2 . l2 FC\ 2 )R2*\ 2 
= -MR 2*1 2 

and (algebra omitted): 

M = Ml M2 

Ml = [V? cia Ral a ] = 1 - ^)KC1) 

M2 = — [Vfl 2 . l2 FC1 2 ] 

= [KCj + 4# [I- KCj ]] 

where for j = 1 , 2, . . . , m 

KCj = kcj 2 kcj 2 
Atj = i j — ioj 
f 0 Ij < toj 

6j = \ 

l 1 £j > ioj 

Here, KCj , I — KCj are projections on and perpen- 
dicular to the cable directions, kcj. If cable Cl is 
inelastic, then Afl = 0 in Ml. 

The kinematic relations u(t>), v(u) can be assembled 
by using equation (116) and then the simulation equa- 
tions can be given by applying equation (9a). Fur- 
ther description of the resulting simulation equations 
is omitted because the treatment here is outside the 
pattern of interest for the applications work; that is, 
there is no subset of the coordinates u that represents 
the inelastic suspension. 

3. The EOMs for v given in reference 16 for the 
two systems (a) and (b) of figure 4 with inelastic sus- 
pensions can be obtained by the procedure outlined 
in section 3, equations (26)-(30), by using generalized 
coordinates which are reference point coordinates like 
those in part a of tables 9 and, part a of table 10, 
except that the reference point is moved to the load 
attachment point for system (a), and to the cable mid- 
point for system (b). 
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Table 9. Simulation equation summary: multicable suspension 


(a) Configuration kinematics: 
f V1 *N\ (I 0 


V — 


VT 


N 


~ Au — 


uj 1 1 

V w2 2 / 

^22 = — S(Rl*a x ) 
^23 = —TV, 2 5 r (i?a 2 * 2 ) 


0 \ /™V\ 


/ A 22 A; 


23 


0 

\o 


I 

0 


0 

/ 


TV, 2 

o 

o } 


wlx 

V &*25 ) 


u ■ 


( F1 V \ 

ulx 
w22 

V ^a2; ) 


( \ 

AI1 t 

V A- J 


V = 


I 

0 

0 


0 

0 

0 


0 

I 

0 


0 

0 

I 


~T 2 , n A 22 S(Ra2* 2 )/ 


\ —T 2 t N T 2iN 

(b) External forces and inertia coupling terms: 
fo = fg + fa - X - D A u 

( EOIat \ /ml g N + FAl jV 

F02iv m2 gw + F A2n — m2 ( A 22 tcdx A 2 3 uj2 2 + TV, 2 F&2* 2 ) 

MOli MAli — F^lx) Jl u)l\ 

\ M02 2 J V MA2 2 - 5(w2 2 ) J2 w 2 2 



f V1 *N \ 


vr N 


Uj\\ 


u2 2 / 


fo = 


A 22 wli = -T Nj1 S(uli) S(Rl * ai ) wli = T/v,i 5 2 (wli) i?l*ai 

A 23 w 2 2 = -T Ni2 [5(w2 2 ) S(Ra2* 2 ) + S(£a2* 2 )]w2 2 = T w , 2 S(u>2 2 ) (S(a;2 2 ) Ra2\ + Ra2* 2 ) 
Tn ,2 Ra,2*2 = T^t 2 5(w2 2 ) Ra2* 2 
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Table 9. Concluded. 


(c) Suspension forces: 



(FC 1 N \ 


( 1 \ 

fc - 

FC2 n 

= A! FC1 n = 

-I 
A 22 

MCh 



M (72 2 ) 


Uia/ 


FC 1 


N 


Inelastic cables: 

FCIn = -[A ,t D- 1 A']- 1 A ,r D- 1 /o 

A ,T D~ l f o — “y •F’OIat - “2 -F02jv + A 22 Jl -1 MOli + A 23 J2" 1 MQ2 2 
A ,t D~ l A' = I + A 22 Jl- 1 A\ 2 + A 23 J2- 1 A 23 

Elastic cables: 


7 =1 ma.i{0, Kj (1 - lojjtj ) + c, |*} Raj 2 
Raj 2 = R2*j2 + ^a2* 2 


FC1 JV = T n<2 £”1 1 
Raj 2 = R2* j 2 + Aa' 
= Raj 2 Ra2* 2 /Pj 


(d) Simulation equations: 

(SF l N \ 

SF 2 n 
SM li 


sf 


D~ l { fo+fc) 


\SM2 2 ) 
u = A~ l sf 


vr N = SF 1 N 


( (FOIat + FCl N )/ml \ 
(F02 N + FC2 N )/m2 
Jl " 1 (MOli +MCI1 
V J2 -1 (MQ2 2 + MC2 2 ) } 


Cc 1 1 = S Af 1 1 
6j2 2 — SM2 2 

Ra2*2 = T 2>n [SF2 n - SFl N - A 22 SMh - A 23 SM2 2 ] 
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Table 10. Simulation equation summary: single-cable suspension 



(a) q c = pitch rotation (about jN) (b) f c = roll rotation (about ic) 

Inertial cable angles. 


T c ,n — Ei{4> c ) Ei(9 c ) 

Vale = {£ e c cos tc, -t 4>C , £) T 



A22 = — TV.i S(Rl*ai) 
A 23 = TV, 2 S(R 2 * 1 2 ) 



B43 = — T c , 7 V A 22 
B44 = — T c ,;v j 4 2 3 
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Table 10. Continued. 


(b) Configuration kinematics for cable angles relative to T\\ 




(a) Pitch rotation (about jl) (b) Roll rotation (about ic) 

Cable angles relative to 7\ • 


T C) i = Ei(A<f> c ) E2{AQ c ) 

AVal c = {£ A 6 C cos A0 C , —£A(j> c , £) T 
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V 2 \ 


I A 22 A 23 Tm,c 
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0/00 



UJ2! ^ 


(00 / 0 / 
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-b Ra\\ — -j- £ kc\ 

A 22 = — Tn } i S(R1*1 1 ) 

A03 = Tjv } 2 S(i? 2 *l 2 ) 
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I 0 0 0 \ 


/ yi V\ 

uli 


0 0/0 


^2V 

uj ( 2.2 


0 0 0 / 


(jll 

AFal c / 


\ ~T c ,n T CtN 5 43 -844 / 


V ^2 2 / 


B43 = ~T CiN A22 
B 44 = — 5c,/v ^23 
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Table 10. Continued. 


(c) Configuration kinematics for cable angles relative to Th' 


kN, kh 


ic' 



iN 


ih 



(a) Heading axes 

Th,N = £ 3 (^ 1 ) 

T c> h — E^ch) E 2 (6 ch ) 

T c ,n = T c>h T hiN 

V alh c = 6 ch cos <f> ch , -l <j> ch , t) T 


(b) Pitch rotation (about jh) 
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N 
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\ValhJ 


(V i* N \ // 0 0 
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ZU)N = 'vtv S(kc N ) kN N - sin <p C h, — cos <p C h sin 0 c h, O) 2 
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( I 000 ^ (V l* N \ 

0 0 7 0 

7 


( viv \ 

a;li 
uj2 2 

\ValhJ 

Baz ~ —T c ,n A 22 
Baa — —T c> n ^23 




0 


\-T c , 


0 
0 
0 

n T CiN 


0 

I 

0 


V? 


N 


a;li 


B43 B44 / \ uj2 2 / 
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(c) Roll rotation (about ic) 



Table 10. Continued. 


(d) External forces and inertia coupling for cable angles relative to T} x \ 
f° = fg + fa — X — D A u 


/foi n \ 


/ ml g N + FAl N \ 

F02 n 


tti2 gw + FA2jy — tyi2 {A 22 wl\ V A 23 tcj2 2 T;y )C Valhc) 

M Oli 


AfAli — 5(wli) Jl a;li 

M02 2 / 


\ M A 2 2 — <f>(t»;2 2 ) J2 uj 2 2 / 


Auxiliary expressions for A: 

/ @ch \ 


0ch 

\iJ 


— [diag(£ cos 0 c /j,, 1)] 1 Valh c 


ioch c — (0c/n ®ch cos <f>chi $ch sin^^) 

'M { 1 tan 9\ sin 0i tan#i cos0i^ 


6 


0 cos 0i 


— sin 0i 


1 


\ipiJ \0 sin 0i /cos 6 1 cos0i/cos#i/ 


uc c = ip i kN c -h u)ch c 

zw N = (t/t + 01 tan 0i) ZU! N + 


kc A y(ujc • kN) — T^ jC ujc c (kc • kN) 


Terms from A u: 

A 22 ^’li == T n>1 S^{lo li) Rl*ci\ — ZiOpf cos$i 0i — ZOJn 0i 6\ 
A 23 ^ 2 2 — ~T n>2 S^{uj 22 ) J? 2 *l 2 
T NiC Valh c = T NiC S(u;c c ) Val h c 



Table 10. Concluded. 


(e) Suspension forces: cable angles relative to Th- 


fc = 


(FC1 n \ 


kCN \ 

FC2 n 

= H TC = 

—ken 

MCh 


€ii 

\MC2 2 ) 


V ~&2 ) 


TC 


£li = S(ftl*ai) Ti t n fccjv 
£2 2 = S(R2*l 2 ) T 2 , n kc N 


( max{ 0, K (£ — t o ) + c £} 


TC = { 


{ ~(H T D~ l H y 1 H T D~ l fo 


(elastic cables) 


(inelastic cables) 


H t D~ l fo = kcj, (F01 n / ml - F02 N /m2) + £ll J l" 1 MOU - J2~ l M02 2 
H t D - 1 H = Jl- 1 *li +&I J r 2" 1 £2 2 


(f) Simulation equations for cable angles relative to Th- 


(SF 1 N \ 


sf = 


SF2 n 

SMh 


= D~ l ( fo + fc ) 


\ SM2 2 J 


ii = A 1 sf 

vr N = sfi n 

wli = SMh 
Cj2 2 = STvI2 2 

Valh c = T CiN [SF2 n - SFl N - A 22 SMh - ^ 23 SM2 2 ] 
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APPENDIX C 

SIMULATION EQUATIONS FOR BIFILAR SUSPENSIONS 


INTRODUCTION 

Several multipoint suspensions were developed or 
considered in the 1965-1975 period of research for 
the Heavy Lift Helicopter, including the inverted-Y 
suspensions (fig. 5(a)) with or without spreader bar 
(ref. 42), or with active longitudinal and lateral con- 
trol of the attachment points by means of control arms 
(ref. 32); the inverted- V suspension (fig. 5(b)) (refs. 15 
and 43), and the three-point suspension (fig. 5(c)) with 
active vertical winching of all three cables and ac- 
tive lateral movement of the forward attachment point 
(ref. 33). The bifilar supension (fig. 5(d)) is of ana- 
lytical interest as a tractable approximation of some 
practical suspensions, and is the case treated in this 
appendix. These suspension designs, along with wind- 
tunnel and flight-test results and additional bibliogra- 
phy, are discussed in the references cited above. 

The object of these multipoint suspensions was to 
stabilize difficult loads developing significant aerody- 
namic specific forces and moments, such as the 8- by 
8- by 20-ft standard cargo container (MILVAN). In 
a single-point suspension, such elongated loads ori- 
ent themselves broadside to the flight path in a max- 
imum drag attitude and become unstable at speeds 
(40-60 knots) that are well below the power- limited 
speed of the helicopter. In addition, the light natural 
damping factor of the load pendulum motion (less than 
0.1) interferes with rapid, precision load placement and 
can result in pilot-induced oscillations at cruise speeds, 
especially in IFR operations. The two-point suspen- 
sions (figs, 5(a) and 5(b)) provide yaw restraint and re- 
strict pitch attitude to achieve stable flight at higher 
speeds in a minimum drag orientation, and the ad- 
dition of active control of appropriate suspension pa- 
rameters increases pendulum damping sufficiently to 
achieve precision load placement and stable flight over 
the helicopter’s power-limited speed range under IFR 
conditions (refs. 32 and 43). 

For many of these suspensions, the equations of mo- 
tion arc difficult to derive. The bifilar suspension is 
tractable by the present application method in which 
coordinates are selected to represent both elastic and 


inelastic suspensions, but the remaining suspensions 
in figure 5 are not. 

A simulation of the inverted-V suspension with elas- 
tic cables, obtained from equations for general elas- 
tic suspensions in which all cables connect the two 
rigid bodies, is reported in reference 15. The behav- 
ior of elastic and inelastic inverted-V suspensions is 
compared in reference 11, where it is noted that the 
variation of the suspension’s yaw-restoring torque with 
relative yawing of the load depends on cable elastic- 
ity, and that this differs significantly between inelastic 
cables and the cables with elastic properties found in 
practice. Thus, cable elasticity must be considered 
in any accurate simulation of yawing motion. Fur- 
ther, one or two cables will collapse during small yaw- 
ing motions of an inelastic inverted-V suspension, and 
this appears intractably complex to simulate as well as 
unrealistic. Approximate equations are given in refer- 
ences 11 and 12, in which the load-suspension motion 
is represented by three angle coordinates (load yaw 
and longitudinal and lateral swing angles) and cable 
stretching is neglected. 

The inverted-Y suspension imposes one constraint 
when inelastic. If the spreader bar is removed, then 
c = 2. No cable or link in the suspension connects 
two rigid bodies, so that force balance at the bar end- 
points must be utilized in the simulation to locate the 
load relative to the helicopter regardless of elasticity. 
This renders the problem analytically difficult, as pre- 
viously seen for the simpler system c in appendix B. 
No simulation equations were found in the open liter- 
ature for this suspension. However, if the spreader bar 
is sufficiently close to the helicopter, then only small 
motions of the bar relative to the helicopter can occur, 
and the system can be approximated as an inverted-V 
suspension. Similarly, if the bar is sufficiently close to 
the load, then only small motions of the bar relative 
to the load can occur and the system can be approxi- 
mated as a bifilar suspension, as is done in references 9 
and 36. 

Simulation equations for the bifilar suspension are 
given in references 13 and 15 for elastic suspensions 
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(a) Inverted Y 
suspension 



(b) Inverted V 
suspension 



(c) Controllable (d ) Bifilar 

3-cable suspension 

suspension 


Parameters: 

ml, Jl,Rl*l l9 Rl*2, 

m2, J2 R2*3 2 , R2*4 2 
JPoI, fo2 


Figure 5. Multipoint suspensions. 


in terms of the rigid-body coordinates, and in refer- 
ences 36 and 44 for inelastic suspensions. Equations 
of motion axe derived in reference 36 from Kane’s equa- 
tions, but the load-suspension motion coordinates ap- 
pear to be ill-conditioned for ordinary small-angle mo- 
tion, as discussed below. Equations for a fixed-base bi- 
filar pendulum are derived in reference 44 from Kane’s 
equations; this pendulum is equivalent to the load- 
suspension portion of the bifilar slung-load system. In 
addition, linearized EOMs are derived in reference 9. 

In this appendix, the methods of this report are ap- 
plied to the bifilar suspension of figure 5(d), by us- 
ing relative motion coordinates for the load-suspension 
subsystem. The selection of generalized coordinates 
is an issue of interest. The inelastic load-suspension 
has four DOFs, which might be chosen variously from 
among the seven cable and load Euler angles, but 
many of these choices are ill-conditioned; that is, the 
relation of dependent motion variables to some sets 
of four DOFs becomes undefined or weakly depen- 
dent during ordinary small-angle motion. The load- 
suspension geometry is studied to determine well- 
conditioned DOFs, and simulation equations are given 
for one such set, including both elastic and inelastic ca- 
bles. The coordinates used here are similar to those 
used in reference 44. An alternative, well-conditioned 
set is used in the linearized equations of reference 8. 


A third set used in the exact equations in reference 36 
is found to be ill-conditioned. 

System and Constraints 

The bifilar suspension is shown in figure 5(d) along 
with the enumeration of the bodies, cables, attachment 
points, and system parameters used in the derivation. 
The load is suspended by two cables, Cl, C2, attached 
at distinct points, R1,R2 on the helicopter, Bl. It is 
assumed that the suspension is uncontrolled (p = 0) 
and that the line segment between load-attachment 
points parallels the load x-axis, R34 = a i2. The 
reference cable lengths, £ol, £o2, refer to the unloaded 
lengths of elastic cables or the fixed lengths of inelastic 
cables. 

If the cables are inelastic, then two independent 
holonomic constraints are imposed on the configura- 
tion position: 



62 


ei(r) = \R2* N + T Nl2 R2*3 2 -Rl* N -T NA Rl*U\ 

= £ol 

(117) 

£2{r) = }R2* N +T Ny2 R2*h-Rl\-T Nil RV 2i| 

= £o2 

so that the inelastic system has 10 degrees of free- 
dom; 

c = 2, d=10 (118) 

Generalized Coordinates of System with Inelas- 
tic Cables 

The generalized velocity coordinates of the inelastic 
system can be taken as the helicopter’s rigid-body ve- 
locities and four coordinates that define the motion of 
the load and suspension: 



The immediate task is to select suitable coordinates, 
ft. These can be chosen from the seven cable and load 
Euler angles, but, as shown next, some possible choices 
of four coordinates are singular or ill-conditioned for 
ordinary small- angle motion of the system. 

The suspension geometry shown in figure 5(d) is 
a nonplanar four-sided figure governed by the vector 
equation 

R12 + R24 - R34 - R13 = 0 (119) 

These vectors can be expressed in terms of their 
lengths and directions as given in table 11. The vec- 
tor directions can be given in terms of the usual cable 
pitch and roll angles and load Euler angles described 
in appendix A. These angles can be taken relative to 
inertial space or helicopter body axes Angles rela- 
tive to T\ are used here for simplicity in representing 
ordinary motion since these have steady values during 
both steady turns and straight-line flight. The angles 
relative to T\ are indicated by appending A to the 
usual notation. That is, 

T c>1 = £i(A0 c ) E 2 (A0 c ) c e {cl, c2) | 

T 2>1 = ExiAfa) E 2 {A0 2 ) E 3 (AtI> 2 ) j 

( 12 °) 

where T Ct i defines cable axes, {ic, jc, kc} such that 
kc is in the cable direction away from the helicopter, 
and ic is in the helicopter vertical plane of (il, kl), as 
sketched in table 11. Equation (119), when expressed 


in T\ coordinates, gives three scalar equations among 
six of the seven cable and load Euler angles in equa- 
tion (120) (load roll angle does not occur), and three 
of these angles can be given in terms of the remaining 
three. However, in simulations these angles are more 
conveniently calculated from their rates. 

The angular rates are related by the time-derivative 
of equation (119) as follows: 

R12 + R24 - R34 - R13 = 0 (121) 

The time-derivatives in equation (121) can be taken in 
any coordinate frame T a according to the general rule 

If R = £ u then a dK/dt = £u + a u b x £n (122) 

where £, u are the length and unit direction vector 
of R; a d()/dt is the time-derivative in T a \ Fb is any 
coordinate frame in which u is a fixed vector; and a u; b 
is the angular velocity of Tb and u relative to F a - In 
this appendix the velocities relative to T\ are indicated 
by appending A to the usual velocity notation: 

AVij = 4 Rij for any line segment, Rij 

dt 

Awa = X (jj a for any reference frame, T a 

Expressions for the velocities relative to T\ are listed 
in part b of table 11, and then equation (121) is ex- 
panded in terms of the load-suspension angular rates. 
To obtain ft, note, first, that the relative load roll rate, 
Ap2, must be one of the coordinates in ft, since it 
does not occur in the velocity equation and is inde- 
pendent of the six load-cable angular rates which do 
occur. Thus, ft cannot be taken as the four cable-angle 
rates. An inspection of the suspension geometry in fig- 
ure 5(d) confirms that the effect of load roll about the 
line segment R34 on the relative position and velocity 
of any point in the load would be undetectable from 
any knowledge of the suspension geometry alone. 

Secondly, scalar equations in only four of the load- 
suspension angular rates can be obtained from dot 
products of equation (121) with any vector that is 
perpendicular to two of the vector coefficients of the 
load-suspension angular rates in equation 121 as ex- 
panded in part b of table 11, that is, from dot prod- 
ucts with kcl, kc2, i2, icl x jc2, icl x j2, 
icl x k2, jcl x k2, jc2 x k2, etc. We consider the 
following relations obtained from dot products, with 
kc2, 12, icl x k2, respectively: 
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Aq2 ~ 


= ( j 2 « 

\k2 • 


kc2 \ /£! cosA0 c i icl#kc2 


kc2 / 


Ar2 + 


a k2 • kc2 


\ . a /^1 jcl.kc2\ • 

J A0cl_ Uk2.kc2 J 


(123) 


A0, 


cl 


■( 


£2 cos A<p C 2 ic2 • i2 
£1 cosA(^i cl icl 


• i2 ) A ^ c2 + ( cos A 4> cl icl • i2 ) A ^ cl ( 


^2 jc2 • i2 


£1 cos Ad (:1 icl 


• i2^ 


A0 c2 (124) 


/£2icl xjc2*k2\ A . / a icl • i2 \ ( 

= (— kcl.k 2 ) A ^ + Ulccl.k2 ) A ’- 2+ ( 


£2 cos A0 C 2 ic2 x icl • k2 
kcl • k2 


A6> c2 (125) 


Use has been made of the triple-scalar and vector- 
product identities to obtain convenient forms. 

A review of the coefficients of the angular velocities 
on the right-hand-side of equation (123) shows that 
these are much smaller than 1. That is, in steady 
flight (hover, straight-line, or turning flight) if kc2 is 
perpendicular to j2, icl, jcl, then Aq2 = 0 and, other- 
wise, for small-angle departures from these conditions, 
then |Ag2| << max{|Ar2|, jA0 c i|, |A0 c i|}, assuming 
£1, a are of similar size. This reflects the suspension’s 
restriction of the load relative pitch. Similarly, using 
the dot product of kcl with equation (121), it can be 
shown that |Ag2| << max{|Ar2|, \A0 c2 \, |A<£ c2 |} near 
the condition that kcl is perpendicular to j2, ic2, jc2. 
Consequently, any use of Aq2 as an independent co- 
ordinate will yield equations for the dependent an- 
gle rates that are singular, or nearly so, in Aq2, so 
that Aq2 cannot be selected as a coordinate of ft. In 
that case, ft cannot contain the load angular velocity, 
Au;2 2 . 

A similar analysis of equation (124) shows that A 6 C \ 
has first-order dependence on A 0 c2 and second-order 
dependence on the cable roll rates. This follows after 
assuming that £1,£2, a are of similar size and that in 
the vicinity of steady flight conditions, i2 is nearly 
perpendicular to jcl, jc2 and that it forms moderate 
angles with icl,ic2. In that case, one of A0 c i, A6 c2 
must be included in ft in order to represent the effects 
of cable pitching motion, but not both, since these are 
nearly mutually dependent coordinates. 

A similar analysis of equation (125) shows that the 
coefficients of A^ c2 , Ar2 are first order and the coef- 
ficient of A# c2 is second order. This follows, assum- 
ing that £l,£2,a are of similar size, that T c \,T c2 are 
nearly parallel frames, and that k2 forms a moderate 


angle with kcl. In that case, any two coordinates 
from {A<^ c i, A <^ c2 , Ar2} can be selected to repre- 
sent the cable roll-load yaw motion. The net result 
is that ft must be selected to contain (1) Ap2, (2) 
two. rates from {A<^ c i, A<£ c2 , Ar2}, and (3) one of 
{A# c i, A0 c2 }; for example, 

ft — (A^d, A<^ c i, Ap2, Ar2) or 

(A^d, A0 c i, A 0 c2 , Ap2) (126) 

The first of these, with load roll and yaw and rear- 
cable angle rates, will be used below for the simulation 
of the relative load motion. Equations for the depen- 
dent rates, A#2, A0 c2 , A0 c2 are included in table 11. 
The load-suspension generalized coordinates selected 
in reference 36 are 

ft = (A0 cl , A0 c i, A $ 2 ) 

where £ is the load pitch angle measured about jl. 
This set includes a load pitch-angle rate and is nearly 
equivalent to the set (A<9 ci , A</> c i, Aq2, Ap2), which 
was shown above to be ill-conditioned. Similarly, the 
set in reference 36 can be shown to be ill-conditioned; 
that is, the equations for load yaw-angle rate and A0 c2 
in terms of ft are singular, or nearly singular, for or- 
dinary small- angle motion. The coordinates used in 
reference 9 are equivalent to the second set in equa- 
tion (126) corresponding to three cable angle rates and 
load roll, and are well- conditioned. The coordinates 
used in reference 44 are equivalent to those used here, 
that is, the direction angle rates of the cable and the 
load roll and yaw rates. The principal difference lies in 
the use of polar coordinates to define cable direction in 
reference 44, where one coordinate is undefined when 
the cable is vertical. 

The angular velocity coordinates considered here do 
not exhaust the possible choices of ft. For example, the 
load-suspension geometry can be viewed as consisting 
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of two noncoplanar triangles with sides (R12, R13) 
and (R24, R34). The orientation of the triangle 
(R12,R13) can be defined by the rear-cable angles, 
the orientation of the triangle (R24, R34) relative to 
the triangle (R12,R13) by the roll angle about the 
common diagonal (R23), and the orientation of the 
load relative to the triangle (R24,R34) by its roll 
angle about (R34). These four coordinates are well- 
conditioned but the resulting i>(u) is complicated. 


Simulation Equations Using Generalized Load- 
Suspension Coordinates Relative to Helicopter 
Body Axes 


The generalized velocity coordinates for the system 
with elastic cables u, are to be chosen as the gener- 
alized coordinates of the inelastic-cable system aug- 
mented by two length rates such that (1) the aug- 
mented set comprises generalized coordinates of the 
system with elastic cables and (2) the length rates be- 
come zero if the cables are inelastic. Here, the two ca- 
ble lengths are the holonomic constraints of the inelas- 
tic cable system, and their length rates can be added 
to (Vn^ciHi) along with one of the well-conditioned 
coordinate sets listed in equation (126) for ft, to obtain 
u. Taking ft to be the rear-cable angle rates and the 
load relative roll and yaw rates, then u can be taken 


as 


( 


vv 


N 


\ 


U)l\ 

AVl3 c i 

h 


(127) 


\ (Ap2, Ar2) r / 

where AV13 is the relative velocity of the rear cable, 
which, using equations from part b of table 11, is given 
in terms of the cable angle and length rates by 


AV13 c i = (n A0 c icosA<ki, -n A 4> cU ll) 1 


The linear velocity components of AV13 C \ are pre- 
ferred as generalized coordinates over the cable angle 
rates A0 c i, A</> c i owing to the simpler kinematic rela- 
tion v(u) obtained. Then the load c.g. inertial velocity 
is given from the derivative of 

R2* = Rl* + Rl*l + R13 + R32* 


as 


V2* = VI* +u;lx Rl*l + (u>l x R13 + AV13) 
— u>2 x R2*3 


or 

V2 * N = Vl\-T Ntl 5(Rl*3i) wl x + T N%c i AVU cl 

+T Ni 2 5(/?2*3 2 ) u;2 2 (128) 

The load angular velocity can be written as 

u;2 2 = T 2j i uj li + Ao;2 2 (129) 

It remains to determine Aq2{u) for the case of elastic 
cables in order to complete the relation v(u). This is 
obtained from the dot product of kc2 with the velocity 
equation for elastic cables; the result is included in 
table 11. 

The configuration kinematics v(u ), u(v) are assem- 
bled in the equation summary, part a of table 12. 
The cable coordinates AV13 c i,£2 in the inverse rela- 
tion, u(v), are obtained by solving equation (128) for 
AV13 c i and solving a similar kinematic relation for 
the forward cable for 12. The result for A” 1 is seen to 
be simpler than for A. Additional equations for the de- 
pendent angular rates Aq2(u) ) Acjc2 c2 (u) and the rel- 
ative angular velocities needed to calculate T 2 j i,T c2> i 
in a simulation are included in the equation sijmmary 
for completeness. 

The external forces and inertia coupling terms fo 
are assembled in part b of table 12. The notation 
F01, . . . ,M02 for the vector elements of fo is intro- 
duced for brevity in later equations. These repre- 
sent the sum of applied forces and moments result- 
ing from weight (mlg,m2g), aerodynamics, and he- 
licopter rotor output (FAl, . . . , MA2), and inertia 
coupling terms from X + A u. The coupling terms 
comprise a large number of scalar terms that are all 
second order in velocity coordinates from u,v , and 
whose computation can be organized in terms of nat- 
ural vectors and matrices as given in the equation 
summary. If expanded to scalar expressions, a large 
number of terms are obtained. This is done in ref- 
erence 36 where the number of such terms exceeds 
300. If the relative motion is assumed sufficiently 
small (|Ao;cl|, |Au>2| << 0.1 rad/sec) and the cables 
are inelastic, then all terms in A u are negligible (of 
the order of 10 -3 rad sec 2 , 10 _3 g) except a term in 
ojI x Rl*2* x cjI from A 22 u;li. This greatly reduces 
the programming and computations required, but is 
inaccurate in representing the dynamics of larger rel- 
ative motions. 
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Equations for the interaction force fc are given in 
part c of table 12. Each cable connects two bodies so 
that the interaction force can be assembled in terms 
of cable directions and tensions as in equation (11) by 

2 

fc = H i T °j = H TC 
0 = 1 

where {Hj} are given in the equation summary. For 
this example, c = m so that H is also a basis of the 
interaction force space and can be used for both elastic 
and inelastic cables. Note that H is identical to —A 


given by the rows 'of A~ l corresponding to £ 1 (third 
component of AFl3 c i) and fa. 

Finally, the simulation equations are listed in paid d 
of table 12. The total specific force s/, due to 
f°+ fc , is assembled, and its vector elements are 
denotedSFl, . . . , SM2 for convenience. The elements 
of u are obtained by expanding A' 1 sf. If the cables 
are inelastic, then £1 and £2 need not be calculated, 
and various terms in A u containing £\,fa are also 
zero in the computation of /o. 


66 



Table 11. Generalized load-suspension velocity coordinates: bifilar suspension 


(a) Suspension geometry: 
R12 + R24 - R13 - R34 = 0 
R13 = J? 1 kcl 
R24 = J? 2 kc2 
R34 =ai2 

Angles relative to T\: 



Suspension geometry. 



(a) Pitch rotation (about jl) (b) Roll rotation (about ic) 


Cable angles relative to Tj . 

T Ct i = E\(A(f> c ) E 2 (A0 C ) c € {cl, c2} 

T 2 , i = Ei(A<t> 2 ) E 2 {A6 2 ) E 3 {Aip 2 ) 

kcl = cos A0d sin A6 c i il - sin A<j> c i jl + cos A0 c i cos A0 c ikl 
kc2 = cos A(f ) c 2 sin A0 c2 il - sin A^ c2 jl + cos A 4> c2 cos A0 c2 kl 
i2 = cos Ai l > 2 cos A d 2 il + sin Aip 2 cos A d 2 jl - sin A 9 2 kl 
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Table 11. Concluded. 


(b) Suspension velocity equation relative to helicopter body axes, T \ : 

AV12 + AV24 - AV13 - AV34 = 0 

Velocities relative to T \ : 

Awcl = A (j) c i icl + A9 cl (cos A0 c i jcl - sin A<t> cl kcl) 

Awc2 = A(j > c 2 ic2 + A 6 C 2 (cos A0 c2 jc2 — sin A^ 2 kc2) 

Ao>2 = Ap2 i2 + Aq2 j2 + Ar2 k2 
AV12 - 0 

AV13 — h kcl + i\ Awcl x kcl 
AV24 = £2 kc2 + £2 Awc2 x kc2 
AV34 = a Aw2 x i2 


Velocity equation for inelastic cables: 

12 {A 8 C 2 cos A^ c2 ic2 - A0 c2 jc2) - £1 (A0 cl cos A<p cl icl - A<j> cl jcl) - a (Ar2 j2 - Aq2 k2) = 0 


/ ]2«kc2 \ 

= k2.kc2 ) Ar2 + rk2.kc2 ( A ^i cos A ^ci icl • kc2 - A0 cl jcl • kc2) 


Suspension angular rates in terms of (A^ c i, A 0 cl , Ar2) for inelastic cables: 

Aq2 

Awe 2 x kc2 = ^(£1 Awcl x kcl + a Aw2 x i2) 

AQ c2 = Awc2 x kc2 • ic2 / cos Ad> c2 
A ^> c 2 = — Awc2 x kc2 • jc2 

Velocity equation and load pitch rate for elastic cables: 

£2 kc2 + £2(A 0 c2 cos A 0 c2 ic2 - A0 c2 jc2) - AV13 - a (Ar2 j2 - Aq2 k2) = 0 

A/t9 — ( i ^kc2 \ Ar o , kc2.AVl3 1 

q ) Ar2 + - k2.kc2 - a~ k2»kc2 ^ 2 


68 


Table 12. Simulation equation summary: bifilar suspension 


(a) Configuration kinematics: 


■ A u — 


'I 0 0 0 0 

I A22 A 2 3 ^24 ^25 

0/000 
0 >142 2I43 A 44 ^45 


AF13 cl 


,(Ap2,Ar2)^ 


AV13 cl = (/I A 6 c i cos A<£ cll -il A^ci, ^1) T 


Auxilliary expressions for A: 
2 O = k2 • kc2 = k 2 j kc 2 i 
z\ — J2 • kc2 = j 2 j kc 2 \ 
zZ = —1 j a zO 
zr — z\ j zO 
z!3 c i — —z£ T c 14 A:c2i 


Submatrices of A: 

A 22 = -T n ,i S(Rl*2\) 

A 23 = TV,i[Ti, c i -f Ti f 2 sj2 zl3l 
A24 = TV.2 sj 2 zZ 
A25 = T Nt 2 [s*2, sfo + zr sj 2 ] 
A 42 = ^2,1 


fll‘2*! = Rl*li + /I Arcli - T\ |2 ^2*3 2 A 43 = j2 2 zl3j x 

si 2 = (R2*3 x i2) 2 = I st column of S(R2*3 2 ) A 44 = j2 2 zZ 
sj 2 = (R2*3 x j2) 2 — 2 nd column of S(R2*3 2 ) ^45 = 

sk 2 = (R2*3 x k2) 2 = 3 rd column of 5(i?2*3 2 ) 


A 43 =10 zr 


TV, 2 — TV,i T\, 2 
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Table 12. Continued. 


( 1 

0 


-1 


0 


V 


u(v) 

! \ 

Ldli 

AV13 cX | = A 

h 

\ (Ap2, Ar2) T } 

Submatrices of .4 _1 : 

T c \ t N = T ci) i Ti } n 

£33 = T cIA S{Rl* h) + S(Rl3 el ) T cU 

^34 = —Td^i Ti t 2 S(R2*3 2 ) 

kc2pj = TjV, i A;c2] 

£12i = (Rl*2 x kc2)i = S(i?l*2i) Jkc2 x 
£22 2 = (R2*4 x kc2) 2 = S{R2*4 2 ) T 2 ,i kc2 x 

1 0 0 


0 0 0 ^ 

0/0 
~T c i,n T c1> . v B 33 B 3 4 

-he 2% kc2jf -£12f £22^ 

o b 53 b 3 4 

/ 


D _ I 1 0 0 1 t-> 

Ss3 - - 1 0 0 1 ) r w 


B54 = 


/ yr "\ 
vt n 
evil 

\ uj2 2 ) 


0 0 1 


Relative angular velocities: 


Ag2 = z£ h + 213^ AV13 c i + zr Ar2 
Ao? 2 2 = (Ap2, Ag2, Ar2) T 

AF24i = (AV13 + AV34)i - T ijCl Ayi3 cl + a (Ar2 j2 x - Aq2 k2 x ) 
A0 c2 = ic2j AV24i / £2 cos A<p c2 
A$ c2 = — jc2j AV^li / ^2 


Au;c2 c2 = (A<^ C 2 , A0 c2 cosA<p C 2 i — A0 C 2 sinA<£ c2 ) T 


Aa;c2i — A0 c2 ic2i + A0 c2 jli = (A$ c2 cos A# c2 , A0 c2 , — A0 c2 sin A0 c2 ) T 
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Table 12. Continued. 


(b) External forces and inertia coupling terms: 


fo = 


fg + fa- 

X - 

/FOl/v \ 

( 

F 02 N 


M Oli 


V M02 2 ) 

V 


m2 gx 4- F A2n — m 2 {A 22 wli + A 23 AkT3 c i + - • _ 1 (2 + A 25 {Ap 2 , Ar2) 
MAli — 5(wli) J1 wli 

MA2o — S(u>2 2 ) J2 uj2 2 — J2 (A42 wli + A43 A V 1 3 c i + A44 (2 

+A 25 (Ap2,Ar2) T )/ 


Auxiliary expressions for A: 
to 2 2 — -^ 2,1 1 “i - Alo2j2 

ooc\\ — col\ + Au;cli 

lo2c22 = Ao;22 — T 2i \ Au?c2i 
u;c 2 cl cl = T c i,i (Aud2 x - Audi) 
iO = kc2 • ( Acz>2 — Au>c2) x k2 
= kc2j T h2 S(to2c2 2 ) k2 2 
zl = kc2 • ( Au;2 - Au;c2) x j2 
= kc2j T h2 S{lo 2 c 2 2 ) j2 2 


Submatrices of A: 

A 22 = -T Ni1 [5(w h) S(R1*2\) + S{R1*2\)} 

A 23 = ^iv,i[5'(wcli) T liC i 

+ T 1>2 (S(w2 2 ) sj 2 zU'A + sj 2 il3&)] 
A 24 = T N , 2 [zi I + z£ S(oj2 2 )]sj 2 

A 2 o = Tn , 2 [/'I2, £*22] 

j[xl2 = S(a;22) SZ2 

\i2 2 = 5(o;2 2 ) (sfc 2 + zr sj 2 ) + SJ2 
A 42 = — 5(Act?22) 2 * 2,1 
A 43 = j 2 2 il3^ 

A 44 = zi j 2 2 


'0 0 
A 45 = 1 0 ir 
0 0 


zt = — z£ zO/zO 
ir = (zl - zr z0)/z0 
il3 cl = -(H kc2 + z£(Awc2 - Aucl) xkc2) c i = [zlfzl I + S(u;c2cl c i)] 2l3 c i 
Al*2*! = T 1>cl AV13 C 1 +Ii,2 S{R2*3 2 ) Aw2 2 
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Table 12. Continued. 


(c) Suspension forces: 


fc=[H 1, H 2] 



fccljy = Tjf'i kc\\ 


f kcltf £c2jv ^ 

— fccl/v —kc2fi/ 
fill fl2i 


'TCI’ 
TC 2 


^-£21 2 — f22 2 j 


fec2jv = Tv, i A;c2j 

fill = (Rl*l x kcl) ! = 5(i?l*lj) kcl 1 


f21 2 = (R2*3 x kcl) 2 = S{R2*3 2 ) T 2 , i fcch 
fl2 x = (Rl*2 x kc2)! = 5(Tr2i) fcc2 x 
f22 2 = (R2*4 x kc2) 2 = S{R2*4 2 ) T 2>1 


Cable tensions for elastic cables: 

TCj = max{0, Kj (ij - f'oj)} j = 1,2 

Cable tensions for inelastic cables: 

( TC2 ) = ~{ HT D- 1 #}- 1 H t D-i fo 

H1 t D~ l fo = kcl^ (F01 N /ml - F02 N /m2) + £ll{ Jl" 1 M01 1 -£21^ J2~ l M02 2 

H2 t D~ l fo = kc2^ (F01 N /ml - F02 N /m2) + fl2f Jl~ l M0l 1 - f22£ J2~ l M02 2 

H\ T D- 1 HI = fil2 + fllf Jl" 1 fllj +f21^ J2- 1 f21 2 

H2 T D~ l H2 = n\2 + fl2f Jl" 1 fl2i +f22^ J2 _1 f22 2 

H1 T D~ l H2 = /i!2 &cl£. fcc2 w + fllf Tl -1 fl2i +f2l|’ J2 _1 f22 2 
where 

/i!2 — (ml + m2)/ml m2 
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Table 12. Concluded. 


(d) Simulation equations: 


/ SF In\ 


sf = 


SF2 n 
SMI i 


D- 1 ( fo+fc ) 


V SM2 2 / 


ii = A 1 sf 
VV N = SF1 n 


ujIi = SM lj 


AV13 cl = T cliJV {SF 2 n - SF1 n ) + B 33 SM1 X + S 3 4 SM 2 2 
12 = -H2 t sf 


(Ap2\ f 1 0 0\ 

\^Ar2 y V° 0 V 


(5M2 2 -T 2 ,! 5M1i) 
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APPENDIX D 

SIMULATION EQUATIONS FOR DUAL-LIFT SYSTEMS 


INTRODUCTION 

Simulation equations are derived in this appendix 
for the three dual-lift configurations shown in figure 6. 

Rigid loads that are sufficiently long to allow safe 
separation of the two helicopters can be suspended di- 
rectly below the helicopters (fig. 6(a)). Several isolated 
civil operations have used this configuration (ref. 4), 
including a load carried with a separation of 1.5 ro- 
tor diameters. Shorter loads can be carried by various 
configurations considered in the early studies (refs. 1 
and 2). Flight tests and control-system analyses have 
focused principally on the arrangement shown in fig- 
ure 6(c), which uses a spreader bar (refs. 3-7, 17, 20, 
and 21). Limited flight tests were conducted using a 
relatively heavy spreader bar that was 2 rotor diame- 
ters long in order to ensure safe helicopter separation. 
Control automation is expected to result in minimum 
separations of 1.25 rotor diameters and a spreader-bar 
weight penalty of about 5% of payload (ref. 4). Al- 
though this weight is small relative to the other masses 
in this system, it is included in the equations of mo- 
tion, and the configuration is represented as a system 
of four rigid bodies. An alternative three-body con- 
figuration without a spreader bar (fig. 6(b)) has also 
been considered in references 18 and 19. The stabiliza- 
tion of difficult loads has not been considered in the 
suspension designs of figure 6 and may lead to new 
suspensions in the future. 

Until recently, work on the equations of motion of 
dual-lift systems was limited to approximate models 
tractable for control studies of these complex systems. 
The cables have been assumed inelastic in all cases. 
In reference 17, the slung- load systems are approxi- 
mated as point masses linked by fixed-length cables 
in tension, and general equations of motion for these 
systems are given from d’Alembert’s principle. The 
general equations required inversion of a d x d ma- 
trix analogous to that in equation (14), but analyt- 
ical inversion did not appear feasible for a 12-DOF 
point-mass model of the system in figure 6(c). In later 
work on this system at Princeton University (refs. 6 
and 7) Lagrange’s equations were used, assuming a 


point-mass load and thin-rod spreader bar (16 DOFs). 
Linearized equations specialized to hover are obtained 
(ref. 7) and the results include a real-time simulation 
of two-dimensional motion (seven DOFs, ref. 6). More 
recently, nonlinear equations of two-dimensional mo- 
tion have been given for a point- mass approximation 
of the three-body system shown in figure 6(b), which 
has four DOFs (ref. 18). General nonlinear EOMs for 
the dual-lift system with spreader bar similar to those 
given here were initially reported in reference 8; non- 
linear EOMs in terms of the rigid body velocity coor- 
dinates for the same system with inelastic suspension 
are given in references 19 and 20. 

Simulation equations are derived in this appendix 
for the three dual-lift configurations of figure 6 using 
the methods of this report. The results accommodate 
elastic or inelastic cables and dissimilar helicopters. 
All three dual-lift systems can be integrated in to a 
single set of simulation equations. In addition, equa- 
tions are given for an approximate model of the dual- 
lift system with spreader bar, using point- mass heli- 
copters and load, and a thin-rod spreader bar. With 
elastic cables, this system has 14 DOFs; with inelastic 
cables it has 10 DOFs. Cables, rigid bodies, attach- 
ment points, and system parameters are enumerated 
in figure 6. The system parameters are the masses, in- 
ertia matrices, and body-axis coordinates of the cable 
attachment points for each body, and the reference ca- 
ble lengths, which refer to the unloaded lengths 

of elastic cables or the fixed lengths of inelastic cables. 

Dual-Lift System for Long Loads 

System and constraints- The configuration in 
figure 6(a) accommodates long loads and consists of 
three rigid bodies Bl, B2, B3 connected by two ca- 
bles, Cl, C2, attached at Rl, . . . , R4. If the cables are 
inelastic then each imposes one holonomic constraint 
on the configuration motion by fixing the distances 
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Parameters: 

ml, J1.R1M, 
m2, J2, R2*2 2 
m3, J3, R3*3 3 , R3*4 3 
jPol, /Po2 


(a) Long loads 



Parameters: 

ml, Jl, Ri*i, 
m2, J2, R2*2 2 
m3, J3, R3*3 3 
/Pol, Ro2 


Parameters: 

ml, Jl, Ri*i, 

m2, J2, R2*2 2 

m3, J3, R3*3 3 , R3*4 3 

m4, J4, R4*54 

/Pol, /Po2, /Po3, ,Po4 = /Po3 


(b) Pendant suspension 


(c) Spread bar suspension 


Figure 6. Dual-lift suspensions and system parameters. 
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R13, R24: 

£l(r) = \R3* N + TV , 3 R3*3 3 - Rl* N 
-T na ill'll | = tol 

> 

£2{r) = \R3* N + T Nt 3 i?3*4 3 - R2* N 
-T n , 2 R2*2 2 \ = £o2 

and then, for inelastic cables; 

c = 2, d = 16 


(130) 


(131) 


Generalized velocity coordinates and configu- 
ration kinematics- The generalized velocity coordi- 
nates u are to be selected so as to separate motion 
caused by cable stretching from motion with invariant 
cable lengths. This is readily done by extension of the 
single-cable example of appendix B. First, assume that 
the c.g. velocity of Bl and the angular velocities of all 
bodies {V1* N ,ul i,u>2 2) w3 3 ) can be included in u, and 
then observe that the c.g. velocities of B2, B3 can be 
given from the derivatives of 


R3* N = R1* N +T Nt i Rl*u + R13 N 

—Tv , 3 -R3*3 3 

R2* n - R3* N + Tv, 3 # 3*4 3 - R2A n 
—Tv , 2 R2*2 2 


(132) 


as 

V3* N = VV N - Tv, i S{R1* li) u)\\ 

TTv,ci V"13 c i + Tv , 3 S(R3*3 3 ) oj3 3 


These velocity coordinates separate the cable rota- 
tional and stretching motion, and, together with 
(yV N ,uli,u2 2 ,w3 3 ), make up 24 generalized coordi- 
nates for the system with elastic cables which contains 
a subset of 16 coordinates for the system with inelastic 
cables. 

The configuration kinematics v(u) are assembled in 
the equation summary of table 13. The 18 x 18 ma- 
trix, A, contains 6 rows representing equations (133), 
and 12 rows from the unit matrix for the coordinates 
that are in both u and v. The reverse relation, i ;(u), is 
readily given after obtaining V13 c \(v),V2A c2 {v) from 
equations (133). Only transformations and skew- 
symmetric matrices representing Coriolis terms occur 
in A, A" 1 . 

External forces and inertia coupling- The ex- 
ternal forces on the configuration and the coupling 
terms, fo, are assembled in part b of table 13. The 
vector elements of fo are denoted (F01, . . . , M03) for 
brevity in later equations. The fo contains the forces 
and moments on the rigid bodies due to weight, aero- 
dynamics, and rotor output, along with the coupling 
terms, X+DA u, which are all second order in velocity 
coordinates from u, v. The submatrices of A require 
only time-derivatives of transformations (appendix A), 
and A u yields terms which are recognized as centrifu- 
gal and Coriolis accelerations of the forms wxwxR 
and a; x V. 


V2* n = V3* n - Tv, 3 S(i?3*4 3 ) u,3 3 

— Tv,c 2 ^24 c2 + Tv,2 S(R2*2 2 ) lo2 2 ) 

(133) 


where V\3 c \,V2A c2 are the inertial cable veloci- 
ties referred to cable axes T c i, T c2 . As described 
in appendix A, these axes are obtained by defin- 
ing inertial pitch and roll angles for each cable 
(tf>dAi), (^c2 } ^c2 ) ? and then cable axes T c i,Tc 2 are 
obtained from the inertial axes by 

T c> v = E\ (0 C ) E 2 {0 c ) c = cl,c2 (134) 

where kcl,kc2 are along the cables, R13,R24; and 
icl,ic2 are in the inertial vertical plane of (iN,kN). 
For these axes 


u c = 4> c ic + 6 C (cos (j) c jc — sin 4> c kc) c = cl, c2 

(135) 

and then 


Vl3d = t\ kcl cl + (wcl x kcl) c i n 

= (£1 6 cl COS 4>c\ i -a 4>cl , ^1) T 

V24 C 2 = (£2 6 c2 cos <j> c2 , -12 4> c2 , £2) T 


} (136) 


Suspension forces- The resultant forces and mo- 
ments applied by the suspension to the configuration 
of rigid bodies, /c, is given in part c of table 13. Here, 
each cable connects two rigid bodies so that fc can be 
given from equation (11) in the case of elastic cables. 
Since the number of cables and constraints is equal, the 
configuration vectors of the cable tensions {TT1 , TT2} , 
given in table 13, comprise a basis of the suspension 
force space which can be used with equations (18) 
and (19) to give fc for inelastic cables. The results 
are given in table 13. This basis is identical to the for- 
mal basis, —A, given by the rows of A -1 corresponding 
to the stretching-motion coordinates £l(v),£2(v). The 
suspension force parameters for inelastic cables are the 
two cable tensions, and the equations for elastic and 
inelastic cables differ only in the calculation of these 
tensions. 

Simulation equations- Finally, the equations for 
u, are assembled in part d of table 13 from A” 1 s/, 
where sf refers to the applied specific forces and mo- 
ments on the configuration due to fo + fc. The vector 
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elements of sf are denoted SF1, . . . , SM3 for conve- 
nience. The equations for VV N) a;l 1? u>2 2 , o;3 3 are 
identical to those usually obtained from the Newton- 
Euler equations. The equations for the (relative) cable 
accelerations V13 c i, V24 c2 , depend solely on differ- 
ences of specific forces and moments on the connected 
bodies. For inelastic cables, the third components of 
yi3 cl ,y24 c2l namely, £l,£2 are theoretically zero. 

Dual-Lift System with Pendant Suspension 

The system shown in figure 6(b) consists of 
three rigid bodies with the load suspended by two ca- 
bles Cl, C 2 connected at a common point R3. The 
load and the hardware attaching it at R3 are regarded 
as a single rigid body, B3. This system is simply a spe- 
cialization of the previous configuration to one with a 
single attachment point at the load, R3. Analysis and 
simulation equations are identical to those of the previ- 
ous case except for replacing R4 with R3 throughout. 
An equation summary is therefore omitted. Equations 
for this system are also included in the results for the 
multilift system with m-helicopters and pendant sus- 
pension that is treated in appendix E. 

Dual-Lift System with Spreader Bar 

System description and constraints- The dual- 
lift system shown in figure 6(c) consists of four rigid 
bodies: the two helicopters, B 1, B2, are connected by 
tether cables, Cl, C2 to a spreader bar B 3, and the 
load is suspended from the spreader bar by the bridle 
cables, C3, C4, attached at a common point, R5. The 
load and the hardware attaching it to R5 axe regarded 
as a single rigid body, B4. 

It is unnecessary to make any specializing assump- 
tions about the system parameters, but it is noted that 
existing designs are characterized by bridle cables of 
equal length with an angle from the spreader bar in 
the range 45° — 60° , with a spreader bar length in the 
range of 1.25 to 1.5 rotor diameters, and a spreader- 
bar weight that is about 5% of the system payload 
(ref. 4). Identical helicopters with equal loading are 
usually considered, but the present results accommo- 
date dissimilar helicopters and unequal loading. 

If the four cables axe inelastic, they impose four in- 
dependent holonomic constraints on the configuration 
by fixing the following distances: 


£l(r) = \R3* N +T NA R3*3 3 

-R1* n -T N}1 RV1 1| = £ol 

£2(r) = \R3* N + T/v )3 R3*4 3 

-R2* N -T Ni2 R2* 2 2 | = £o2 

(137) 

£3(r) = +T/v >4 -R4*5 4 

— i?3^ — T/V3 i?3*33| = £o3 

£4 (r) = \R4* n +T na R4* 5 4 

—R3* N — Tjv,3 jR3*4 3 | = £o4 

Then, for inelastic cables, 

c = 4, d = 20 (138) 

Generalized velocity coordinates and configu- 
ration kinematics- The generalized velocity coor- 
dinates u are to be selected so as to separate mo- 
tion caused by cable stretching (four coordinates) 
from motion caused by invariant cable lengths (20 co- 
ordinates). First, note that the subsystem that 
consists of helicopters, tether cables, and spreader 
bar is identical to the system in figure 6(a) and 
can be represented by the identical 18 coordinates 
(V1* N , ul\, u; 2 2 , ^3 3 , F13 c1 , V24 c2 ). As before, in- 
ertial cable angles 0 C \), (0 C 2? #02) and cable axes 
T cU F C 2 are defined (eq. (134)), and then the iner- 
tial velocities, Vl3 c i, V24 c2 separate the tether cable 
stretching and rotational motions (eq. (136)), and the 
rigid body velocities V3* N ,V2* N are given in terms of 
the generalized coordinates by equation (133). 

Secondly, the load motion remains to be given in 
suitable coordinates. The load c.g, velocity is obtained 
from the derivative of 

R4* n = R3* N + R3*b N - T NA /M*5 4 

as 

V4* n = V3* n + V3*5 n + T na S(R4*5 4 ) ll>4 4 (139) 

Assume that u>4 4 can be included in u. It then remains 
to develop appropriate coordinates for V3*5jv. This 
development is given in table 14. 

The inertial velocity, V3*5, results from the iner- 
tial velocities of the two cables C3,C4. These cables 
always form a triangle with the spreader-bar longitudi- 
nal axis i3, and axes can be attached to the plane of 
the triangle with longitudinal axis along the spreader 
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bar, i3, and normal axis, kt, perpendicular to i3 in 
the plane of <73, C 4. This plane can have any roll 
angle, 4 > t , about i3, and then R3*5 has components 
(£3x, tZz) in this plane. For elastic cables, (£3x, CZz) 
are arbitrary; for inelastic cables, they have fixed val- 
ues. The coordinates {0t, £3x, £3 z) are appropriate ve- 
locity coordinates which separate the stretching and 
rotational motion of R3*5. To isolate these coordi- 
nates in V3*5, define the modified spreader-bar ref- 
erence frame F 3 ', from the spreader-bar Euler angles 
{^ 3 ,^ 3 }, and then give the inertial velocity of R3*5 
as the sum of its velocity relative to F 3 ' (VT in part a 
of table 14) and the effect of the inertial rotation of 
F 3 < (AVT in part a of table 14). As seen in table 14, 
the motion of R3*5 caused by cable stretching and 4> t 
is isolated in VT t and the remaining motion, AVT 
depends only on the spreader-bar angular velocity. 

Appropriate generalized coordinates for the system 
can now be given by augmenting the 18 coordinates 
previously used for the system of figure 6(a) with 
(VT t , o>4 4 ); V 4* n is given in terms of these coordinates 

by 

V4* n = Vr N + T N ,t ( VT t + BT w3 3 ) 

+Tw, 4 S{R4'5 4 ) w4 4 (140) 

Finally, the kinematic relations v(u),u(v) are assem- 
bled in table 14, by using equations (133) and (140) 
to obtain V2* N {u), VZ* N (u), V4* N {u) and then re- 
arranging these same equations to obtain VT3 c i(v), 
V24 c2 (v), VT t (v). The remaining 15 rows of A, A -1 
are all from the unit matrix, since the corresponding 
coordinates are in both u, v. The nontrivial submatri- 
ces of A, A -1 are all coordinate transformations and 
skew-symmetric matrices representing Coriolis terms. 

Applied forces and inertia coupling- These 
terms, /o, are assembled in part c of table 14. The 
vector elements of f o are denoted F01, . . . , M04 for 
brevity in later equations. The configuration vector 
fo contains the sum of applied forces due to weight 
(ml g, ...,m4 g), aerodynamic and helicopter rotor 
forces and moments (FA 1 , . . . , MA4), and the inertia 
coupling X 4~ A u. The submatrices of A are defined 
in part b of table 14 except that BT is given in part a. 
Their time-derivatives are obtained principally from 
derivatives of coordinate transformations. A general 
expression from appendix A is repeated here along 
with expressions for the cable-axis angular velocities 
in the required reference frames, wcl c i,wc 2 C 2 , u>t t , and 
expressions for the nontrivial terms of A u. The time- 
derivative of BT is obtained routinely, but a formula 
is omitted for brevity. 


Suspension forces- Equations for the suspension 
force fc are assembled in part d of table 14. Each 
of the four cables <71, <72, <73, C4 connects two rigid 
bodies so that fc can be given as in equation (11) by 

= H i TC i = HTC ( 141 ) 

j = 1 

The matrix H is given in table 14. Its columns also 
form a basis of the suspension force space since c = 
m = 4, so that equation (141) can be used for both 
elastic and inelastic cables, the cable tensions being 
calculated from cable stretch or from equations (18) 
and (19) in the two cases, respectively. Expressions for 
both cases are fisted in the table. The table includes 
formulas for calculating cable line segments from u, 
and these suffice to determine H and the elastic cable 
tensions. 

Alternatively, a basis can be given as the rows of 
A~^ corresponding to the cable-stretching coordinates 
(rows 6, 9, 10,12). This basis is given as A in part d 
of table 14. The corresponding constraint force pa- 
rameters s can be identified as T<71, TC2 , and the J~t 
components of the suspension force, FT, which acts 
on the triangle at R5. Note that these components 
can be given from the cable tensions and, conversely, 
so that A can also be used to assemble fc for both 
elastic and inelastic cables. 

Simulation equations- The total specific force ap- 
plied to the configuration, sf, is assembled in part e 
of table 14 from the vector elements of the external 
forces fo and the interaction forces fc. The vector 
elements of s f are denoted SF 1 , . . . , SM4 for conve- 
nience. Finally, the simulation equations for the vector 
elements of u are given by expanding u = A -1 sf , and 
these represent either elastic or inelastic cables. For 
inelastic cables, the four cable-stretching coordinates 
(coordinates 6, 9, 10, 12 of u) are all theoretically zero. 

Remarks- The three dual-lift systems can be inte- 
grated into a single simulation. The three-body system 
in figure 6(a) is a subsystem of the four-body system 
c. That is, system a is obtained from system c by 
dropping the load and bridle cables and by regarding 
the spreader bar as the load and assigning it appropri- 
ate parameter values. System a can be represented 
by a subset of the generalized velocity coordinates 
and equations for system c, obtained by deleting the 
six load-triangle coordinates and the load forces and 
moments. Further, as previously noted, the pendant 
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suspension system in figure 6(b) is a simple specializa- 
tion of system a to coincident load- attachment points 

(R4 = R3). 

The EOMs for the rigid-body velocities, given in 
reference 20 for dual lift with spreader bar and inelastic 
suspension can be obtained by the procedure outlined 
in section 3, equations (142)- (146), by using general- 
ized coordinates, u, which are reference point coordi- 
nates like those in part b of table 14, except that the 
reference point is taken at the spreader bar c.g. These 
coordinates also result in a simpler kinematic relation, 
v = An , than in part b of table 14 and are of interest 
as an alternate set of generalized coordinates for the 
simulation equations. 

Dual-Lift System with Spreader Bar: 
Degenerate-Body Approximation 

In this section the system is represented by point- 
mass helicopters and load, and a rigid-rod spreader 
bar, as illustrated in part a of table 15. The system 
with elastic cables has 14 DOFs. Inelastic cables im- 
pose four holonomic constraints as in equation (137), 
in which case d = 10. 

The formulation of EOMs for systems with degen- 
erate bodies is outlined in section 5. To account for 
the point masses, 51,52,54, it is only necessary to 
modify the equations in table 14 to remove the angu- 
lar velocities from v, u, along with the associated rows 
and columns of A, A -1 , and the associated moments 
in fo, fc. A 15-DOF system is obtained. 

The rigid-rod spreader bar has two attitude degrees 
of freedom, represented by its inertial heading and 
pitch angles. The spreader-bar body axes are now 
coincident with the special axis frame Fy, defined in 
table 14, and the attachment points on the spreader 
bar are on the i3-axis (table 15). 

To account for the rigid-rod spreader bar 53, define 
the reduced configuration velocity v, in which w 3 3 is 
replaced by the spreader-bar’s pitch and heading rates, 

a3 = ( 03 , ipz) T , and make the same replacement in the 
generalized coordinates, u : 

v = (Vl* Ny V2* N , V3* n , V4* n , a3) r 

u=(Vl* N ,V13 cU V24 c2) VT t , ^3 ) T 

The required kinematic relation, 


is expanded in table 15. The cable velocities 
V13, V24, VT were defined previously in equa- 
tion (136) and part a of table 14. The 14 coordinates of 
u contain four scalar coordinates that 

represent cable stretching and 10 coordinates that rep- 
resent motion with invariant cable lengths. The kine- 
matic relations v(u),u{v) are obtained by specializ- 
ing equations (133) and (140) to the degenerate bod- 
ies. The results are assembled in part b of table 15. 
The submatrices of A, A are coordinate transforma- 
tions and cross products with the spreader-bar’s axes 
of pitch and heading rotations. 

Equations for the applied forces and inertia reactions 
fo , are given in part c of table 15. The configuration 
vector fo is expanded routinely except to note that A 
in Au is from the relation v(u) = W A u, where W 
maps Euler-angle rates to angular velocities (eq. (32)). 
The inertia reactions consist of Coriolis velocities due 
to transformation rates, centrifugal accelerations due 
to spreader-bar angular velocity, and an effect, zl, of 
the Euler-angle rate coordinates. 

The suspension forces (table 15) can be given in 
terms of cable tensions as in table 14, fc= H TC , 
where H here is obtained from H in table 14 by delet- 
ing the elements that generate moments on the point- 
mass bodies. The moment action vectors {£3j} were 
defined previously in table 14 along with equations 
for calculating cable line segments. For elastic cables, 
TC is calculated from the cable lengths as usual. For 
inelastic cables, first obtain the reduced basis in 14- 
dimensional space, 5 = W T H , and then the cable 
tensions are given by equation (37). Inversion of a 

4x4 matrix ( H D~ l H) is again required. Expres- 
sions for its elements are routinely obtained by expan- 
sion of the matrix. The alternative basis, A from the 
rows of A~ l corresponding to the cable stretching co- 
ordinates, can also be used in equation (37) to obtain 

fc = A s = W fc. The cable tensions can be ob- 
tained from s as given in part d of table 14. 

Finally, the simulation equations (36) are assem- 
bled in part c of table 15. The reduced specific force, 

— x j* 

sf = D W (fo + /c), contains the specific forces, 
SF1, . . . , SF4, and t he reduced specific moment on the 
spreader bar, SM3; the latter contains the compo- 
nents of the moment sum, M03 + MC3, along the 
pitch and heading axes of rotation, j3,kN. Last, 
the simulation equations for ii are given by expand- 
ing equation (36). 


Table 13. Simulation equation summary: dual lift for long loads 


(a) Kinematics: 



K13 c i = (£1 0 C \ cos^d, -i\ d>c i, h) T 
V24 c2 = {(-2. 6 c2 cos (t> c2 , -12 <j> C 2 , £2) T 


A 04 — — I’v, 2 S(Rl*li) ^24 — Tc 1 , ,V ^24 

A 25 = Tjv , 2 S{R2*2 2 ) B 26 = -T c i,at ^36 

4 2 6 = —Tn,3 S(i?34 3 ) B35 = T c2) iv -4 2 5 

A 36 = IV.3 S{R3*3 2 ) B 2 6 — —T c2 js! Tfr\3 S(R3 4 3 ) 
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Table 13. Continued. 


(b) External forces and inertia couoling: 
fo= fa - X - D A u 4- fg 


/ F01 n \ 


F02 n 


fo = 


F 03* 
AfOli 


MO 2 2 
\ M03 3 / 


/mlg N + FAl N \ 

m2 + F A2n — ra2(T/v iC i V 13 c j — T/v ?c 2 V r 24 c2 + Aoa ^li 

+^25 W 2 2 + -4.26 ^3 3 ) 

^i3 /j'iv 4- FA3n — 77i3(T/v,ci V 13 c i 4- ^4 24 wli 4* ^4 3 6 <^3 3 ) 
AfAl\ — «S* 1 1 ) 

MA22 — S(uj2 2 ) J2u?2 2 

V MA3 3 - S(u;3 3 ) J3 u;3 3 / 


^A r ,c = 7!v,c »5(u^C c ) 

u^c c = (<^C 7 cos 0c j — 9 c sin<f) c ) T 


c = cl , c2 


i 24 Wli =Tat (1 S 2 Mi) R1*U 
A20 ^2 2 = ^-2V l 2 5 2 (w22) /?2*2 2 
^36 ^3 3 = — T/v , 3 S 2 (u;3 3 ) #3*3 3 
^26 ^3 3 = Ttv )3 S 2 (u>3 3 ) /£34 3 
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Table 13. Continued. 


(c) Suspension forces: 


/FClivN 


{ kcl i v 

0 

\ 


FC2 n 


0 

fcc2^v 



FC 3* 


— feci jv 

—kc2w 


/TCI 

MCli 

= HI TCI + F2 TC2 = 

fill 

0 


\TC2 

M C2 2 


0 

£22 2 



MCh > 


-€31s 

“£32 3 

y 



£llx = (Rl*l x kcl)i = S(Rl*li) T 1jN kcl N 
£313 = (R3*3 x kcl) 3 = S(i?3*3 3 ) T 3jN kcl N 
£22 2 = (R2*2 x kc2) 2 = 5(/?2*2 2 ) T 2 , N A:c2 N 
£32 3 = ( 113*4 x kc2) 3 = S(i?3*4 3 ) T 3) at /cc2^ 

Elastic cable tensions: 

TCj = max{ 0, Kj (tj - toj)} j = 1,2 

Inelastic cable tensions: 

(rC2) = _(ffT 23-1 ^ /0 

Hf D - 1 /o = fccl£ + S11J Jl- 1 MOU - €3lJ- J3- 1 M03 3 

i/J D- 1 fo = kc2 T N + £22^ J2- 1 M02 a - £32^ J3" 1 Af03 a 

Hi D~ l Hr = + £llf Jl-i Clli + €31?’ J3- 1 £31 3 

Hi D~ l H 2 = + £22j J2~ l £22 2 + £32^ J3" 1 £32 3 

Hi D- 1 H 2 = fccl£ fcc2 JV /m3 + €3l|’ J3" 1 £32 3 
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Table 13. Concluded. 


(d) Simulation equations: 


/ SF1 n \ 


/ (F01 n + FCl A 0/ml \ 

SF 2 n 


(F02^ + FC2 N )/m2 

SF3tf 


(F03yv + FC3^)/m,3 

SMli 


Jl _1 (A/01i + MCli) 

SM 2 2 


J2~ 1 (M02 2 + MC2 2 ) 

V SM3 3 ) 


V J3 _1 (A/03 3 + A/C3 3 ) / 


u = A 1 sf 
V1* N = SF1 

F13 ci = T c i iN (SF3 N — SFl N - A 24 SMli ~ ^4 36 SM3 3 ) 

V24 C 2 = T C 2 f N(SF3iv — SF2w + -4.25 SM22 — T/y 3 S(/?3*4 3 ) SM3$) 
Cj\\ = S A/ 1 [ 
u2 2 — SM2 2 
w3 3 = SM3$ 



Table 14. Simulation equation summary: dual-lift system with spreader bar 


(a) Generalized velocity coordinates for bridle cables (C3,C4): 
T t = {i3, jt, kt} 

i3 = spreader bar longitudinal axis, along R43 
kt = perpendicular to i3 in plane of C3,C4 
723*5* = (73x,0,73z) T 
J- 3 , = {i3J3',k3'} 

= £2 ($ 3 ) £3 
T% t 3 ' = E\{(f)t) 

V3*5 t = Tt,N R3*3 N = VTt + A VT t 

R3*5 3 > —T 3 ',t R3*5 t = {£3x, -I3z sinept, (3 z cos<p t ) T 

VT t = T t , 3 ' R3*5y = (tx, -t3z<j> u £3 z) T 

u>3' = ip 3 kN + 6-i j3' = [j3' j3' + kN k3'/ cos0 3 ] • w3 

AVT = — R3*5 x w3' 

A VTt = DT w3 3 


/° 


BT = 


0 


Vo 


73 2 cosA<j) 73z sinA4> \ 

— 73x sinA<f> + 73 2 $m0 3 2<m0 3 73x cosAcj) -f 73z cos<£ 3 £an# 3 

—£3xcosA<f) — 73x sinAcj) / 


where 



Table 14. Continued. 


(b) Configuration kinematics: 



A 2 S = -T Ntl S(Rl* li) 

A 2 c = T Nf2 S(R2*2 2 ) 

A 37 - T N , 3 S(R3*3 3 ) 

A' 27 = T Nt3 S(R3*4 3 ) 

A 2 7 = ^37 ” ^27 ” ^iV,3 5(/?343) 

A 47 = ^37 + Tn^BT 

A48 = T^ 4 S(R4*5 4 ) 



B 25 = — T c i t jvA25 

#36 = ^ 1 c2 1 A r ^26 
#27 ~ — 

$37 = — 7^2, AT ^27 
#48 = —T^nA 43 
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Table 14. Continued. 



T n ,c = T n , c S{wc c ) c = 1 , 2 , 3, 4 , cl, c2, t 

wc c = ( <}> c , 9 C cos 0 C , -0 C sin <f> c ) T c= cl, c2 

- / 
o>t = i3 (fit ~h ^3 

1 — sinfotanQz — cos^tanO^ 

0 cosAcj) sinAcj) 

0 —sinA<j> cosA(f> 

A 2 5 uli = -T n , 1 5(«li)S(ill*li) wli = IV, i 5 2 (wli) m*l! 

A 26 w22 = —TW , 2 S 2 (u>2 2 )R2*22 

A 27 uj33 = —TV, 3 S 2 (w33)(i?3*33 — /?3*43) 

i 37 w3 3 = -TV, 3 5 2 (w3 3 )i?3*3 3 

A , 7 W 3 3 = -T Ni 3 5 2 (w3 3 ) fl3*3 s + T N , t S(ut t ) BT u>3 3 + TV, t BT u;3 3 

^48 w 44 = —Tv , 4 S 2 (yj4 i )R4*O i 
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Table 14. Continued. 


(d) Suspension force: 


/ kcljv 

0 

0 

° ^ 


0 

kc2w 

0 

0 


—kcl 

—kc 2jv 

A:c3at 

fcc4jv 


/TCI \ 

0 

0 

— kc3rf 

— kc4w 


TC 2 

£iii 

0 

0 

0 


TC3 

0 

£22 2 

0 

0 


\TC4/ 

-£31 3 

—£323 

£33 3 

£34 3 


\ 0 

0 

-£43 4 

-£44 4 j 



/ fccljv 

0 

0 

0 ^ 


0 

kc 2jv 

0 

0 


— fccl/y 

—kc2iy 

*3 AT 

ktw 


/TCI \ 

0 

0 

— ?3jv 

—kt?<r 


TC 2 

£lli 

0 

0 

0 


FTx 

0 

£22 2 

0 

0 


\ FTzJ 

-£3i 3 

-£32 3 

£3x 3 

£3*3 


\ o 

0 

-£4a- 4 

-&Za) 



TC 1, . . . , TC 4 = cable tensions for Cl, . . . , C4 


Ft = 

: kc3 TC3 

+ 

kc4 TC4 

= i3 FTx 4- kt FTz 





ai 

= Rl*l 

X 

kcl 

^33 

= R3*3 

X 

kc3 

£3x 

= £3z jt 


>/>► 

GO 

= R3*3 

X 

kcl 

^43 

= R4*5 

X 

kc3 

£4x 

= R4*5 x 

i3 

422 

= R2*2 

X 

kc2 

£34 

= R3*4 

X 

kc4 

£3z 

= —£3x jt 


£32 

- R3*4 

X 

kc2 

£44 

= R4*5 

X 

kc4 

£4z 

= R4*5 x 

kt 
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Table 14. Continued. 


jR13jv = Cl kc 1 N = / TV, ci Vl3 c i dt 
R24 N = £2 kc2 N = J T N)C2 V24 c2 dt 
R3*5 N = f Tjv.t (VT t + A VT t ) dt 
R35 N = 13 kc3 N = R3*5 n + / T N , 3 S(R3*3 3 ) w3 3 dt 
R45 N = £ 4 kc4 N = R3*5 N + J T Ni3 S(R3*4 3 ) w3 3 dt 
Elastic cable tensions: 

TCj = max{0, Kj (£j - toj)} j = 1, . . . , 4 

Constraint force parameters for inelastic cables: 
s = TC = -\H T D~ l H]~ l H T D~ l fo 

H\ t D~ l fo = kcljj ( F01n/tti1 — F03w/m3) + ^ll 1 T Jl _1 M01x — £31 3 T J3 1 M03 3 

H2 T D~ l f 0 = kc2jj {F02 N /m2 - F03 N /m3) + £22/ J2~ 1 M02 2 - £32fJ3- 1 M03 3 

H3 T D- 1 f 0 = kc3 T N (F03jv/m3 - F04Ar/m4) + £33 3 T J3 _1 M03 3 - £43 4 J4 _1 Af04 4 

H4 T D~ 1 f 0 = kc4‘ N (F03 N /m3 - F04 N /m4) + £34 3 r J3- 1 M03 3 - f44 l i J4~ l MM i 

H^D^Hl = im 3 + £11i T J1~ 1 £11i + £31 3 T «73 -1 £31 3 

H2 t D~ 1 H2 = n 23 + £22 2 T J2~ l Z22 2 + £32j J3~^32 3 

H3 t D~ 1 H3 = ^34 + £33 3 r ,/3 _1 £33 3 + £43 lJ4~^43 4 

H4 t D- 1 H4 = ^34 +£34 3 T J3- 1 £34 3 +£44lJ4“ 1 £44 4 

Hl T D~ 1 H2 = kclj, kc2 N /m3 + £31 3 T J3 _1 £32 3 

H\ T D~ l H3 = ~{kcl T N fcc3 w /m3 + £31 3 T J r 3 -1 £33 3 ) 

Hl T D- l H4 = -{kclj, fcc4 w /m3 + £31 3 T J3- 1 £34 3 ) 

H2 t D~ 1 H3 = -{kc2% kc3 N /m3 + £32 3 T J3 _1 £33 3 ) 

H2 T D- 1 H4 = ~{kc2j, jkc4 J v/m3 + £32 3 T J3- 1 £34 3 ) 

H3 t D~ l H 4 = n 3 4 kc3jj kc4 n + £33 3 T J3 _1 £34 3 + £43jJ4 *£44 4 
where 

yi i 3 = (ml + m3) /ml m3 
^23 = (m2 4- ml)/ m2 m3 
/i 34 = (m3 + m4)/m3 m4 
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Table 14. Concluded. 


(e) Simulation equations: 



V13 c i = T chN [SF3 N - SF 1 N - A 25 SMI i - A 37 SM3 3 ] 
V24 c2 = T c2>n [SF3 n - SF2 n + yl 26 SM2 2 - A' 27 SM3 3 ] 
VT t = T t , N [SF4 N - SF3 n - A 48 SM4 4 ] - BT SM3 3 
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Table 15. Simulation Equations: Degenerate-Body Approximation of 
Dual-Lift System with Spreader Bar 

(a) System and spreader-bar kinematics: 



Degenerate-body dual-lift system. 


Parameters: 


ml 

m2 


m3, J3 = J3' diag {0,1,1} 
m4 

Jfol, fo2, fo3, /o4 = jfo3, R3*3, R3*4 


i3 = spreader-bar longitudinal axis 

(R34,R3*4, R3*3) = (-#34 i3, — 7?3*4 i3,+#3*3 i3) 

T^,n = #2(^3) £3(^3) 

Tt, 3 = 


w3 3 = (j 3 3 , 


kN a ) 



/° 

- sin 0 3 \ 

f0 3 \ 


• = 

1 

0 

\t hJ 

Vo 

COS 63 / 




= ^ 33^3 
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Table 15. Continued. 


(b) Kinematics: 



( > 


/i 

0 

0 

0 

0 


^ 2 V 


1 

TV, cl 

1 

2 

r> 

to 

0 A25 

V = 


= Au — 

1 

TV, ci 

0 

0 A35 


K 4 V 


1 

TV, ci 

0 

7 V,t A45 


^3 / 


lo 

0 

0 

0 I ) 


i W N > 


T 

0 

0 0 


^ 13 el 


-T c \,N 

0 

T c \ y N 0 

u — 

V 24 c2 

= A 1 v - 

= 

0 

-T c 2 ,N 

T C 2 ^N 0 


VT t 


0 

0 

-Ti,n T i}i 


V o 3 / 


\ 0 

0 

0 0 


V 13 C1 

V 24 c2 

VT t 

V H3 ) 

0 

B25 
B35 


( Fl " \ 


V 2 ‘ 


V 3 ‘ 


V4' 


N 


N 


N 


Auxiliary expressions: 

^25 = +T/v, 3 5(283) = ( +A: 3 yv, — cos 63 j 3 at ) 


\ a3 / 

f + sin 63 cos ^3 + cos O3 sin ^3 ' 

-fsin# 3 sin ^3 — cos# 3 cos^ 3 

+ cos 63 0 y 


Submatrices of A, A 1 : 

A 25 - -T n> 3 5(i?34 3 ) W3 3 = i?34 A 25 
A 35 = TVs 5(,R3 + 3 3 ) IF3 3 = R3*3 A' 25 

-T Nt 3 S{-R3*3 3 + i?3*5 3 ) W3 3 


A45 = — T/v i3 S(R 35 3 ) W 3 3 

B 25 = — A 35 

£35 = -r c2 ,^ Tjv, 3 5(/?3*4 3 ) TF3 3 = /?3*4 r c2iN A^ 5 
AVT t = — (R3*5 x u,3')* = o3 

/ ^3z cos$* £3z sintfo cos 63 

BT t — —S(R 3 * 5 t ) T tf 3 W3 3 = — £ 3 x sin (f> t £ 3 z sin0 3 4- £ 3 x cos<£t cos# 3 

\-£ 3 x cos (p t -£ 3 x sin <f> t cos 0 3 
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Table 15. Continued. 


(c) Applied forces and inertia reactions: 


f° = f9 + f a ~ X — D A u 
A=WA = diag{I, /, I , J, W3 3 } H 



/ F01/v \ 




F02 w 


m2 + FA2 at — m2(T/v >cl VT3 cl — T/v , C 2 ^24 C 2 + A 25 <*3) 

fo — 

F03 w 

= 

m3 gw + F^43/v — m3(Tjv,ci ^13 c i + ^35 a3) 


F04 n 


m4 4- FA4n - m4(7V,ci ^13 c i + T /v,t VT t + A 45 a3) 


M03 3 / 


\ MA3 3 - 5(c^3 3 ) J3 u>3 3 - J3#3 3 S3 / 


Transformation rates: 

Tn,c = Tn iC S (u ;c c ) c = cl, c2, t 

wc c = (0 C1 (9 C cos 0c, -0 C sin0c) T 


c = cl , c 2 


/ 


— 


0 


- sin 0 3 \ 


f*‘\ 




W£ t a£ 


0 cos (fit sin (fit cos #3 
\ 0 — sin (fi t cos (fi t cos $3 / v 03 y 

Other submatrices of .4 n: 

2I3 = W 3 3 a 3 = - 0 3 103 (cos 0 3 , 0 , sin 0 3 ) r 
= A25 = T Nf3 [S 2 (u; 33 ) + 5(2:13)] i 3 3 

/— 0 | COS 0 3 COS 03 + 203 0 3 sin 0 3 sin 03“ 03 COS 0 3 COS 03 ^ 
= —0 2 cos 0 3 sin 03 — 203 0 3 sin0 3 cos0 3 — 0| cos0 3 sin0 3 

\ —0 2 sin 0 3 

A 25 o3 = i?34 z2n 
^35 oc 3 = i?3*3 z2 j\ r 

A 4 5 o3 = T w , 3 [( 5 2 (w 33 ) + S(zl 3 )) R 35 3 + S(u>3 3 ) F35 3 ] 
where 

R35 3 = — F3*3 3 + F3*5 3 
R35 3 = T 3 , VT t 


J 


93 

( 2 - 2 - 



Table 15. Continued. 


(d) Suspension forces: 



/ kclpi 

0 

0 

0 \ 


0 

kc2w 

0 

0 

fc = H TC = 

— kc \ tv 

—kc2^ 

kc3w 

A;c4/v 


0 

0 

— fcc3/v 

— fcc4jv 


V -£31 3 

-£323 

£33 3 

£34 3 / 

Elastic cables: 




TCj = max{0 } 

Kj (ij - 

- to 3)} 

J — li - • 

-,4 

Inelastic cables: 





/ A:c1tv 

0 

0 

0 ^ 


0 

kc2jv 

0 

0 

H = W T H ~ 

— fccl /V 

—kc2pj 

fcc3/v 

fcc4;v 


0 

0 

— fcc3/v 

— A:c4;v 


^ — £31 3 

— £32 3 

+.£33 a 

+£34 3 / 

&j 3 j = 




TC = -\fi T D 

_1 H\~ l ' 

77 r D' 1 

W T /o 



/TCI \ 
TC 2 
TC 3 
\TC4/ 


w = /, /, /, TF3 3 } 

D = iy T D W ~ diag{ml 7, m 2 7, m3 7, m 4 7, J 3 } 
J3 = W^3^ J3 TF 33 = J3' cftag{l, cos 2 ^ 3 } 
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Table 15. Concluded. 


(e) Simulation equations: 


sf = D 1 W T ( fo + fc) = 


/ SF ljv\ 


/ (F01 N + FCl N )/ml \ 

SF2 n 


(F02jv + FC2n)/tti2 

SF3 n 

= 

(F03 N + FC3jv)/m3 

SF An 


(F04^ + FC4n)/tti4 

SM3 ) 


J3 _1 W3 3 T (M03 3 + MC3 3 )J 


it = A 1 sf 

VVjf = sfi n 

Vl3 cl = T chN {SF3 N - SFl N - R3'3 A' 25 SM3) 
V24 c 2 = T C 2'N(SF3n — SF2n + R3*4 A SM3) 
VT t = T t , n (SF4 N - SF3 n ) - WT t SM3 
a3 = SM3 
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APPENDIX E 

SIMULATION EQUATIONS FOR A MULTILIFT SYSTEM WITH PENDANT 

SUSPENSION 


INTRODUCTION 

This appendix gives simulation equations for the 
multilift system shown in figure 7, in which m heli- 
copters support a load with a single cable attaching 
the load to each helicopter. This system was sug- 
gested in discussions by P. K. A. Menon (ref. 18) as 
one of potential interest using remotely controlled heli- 
copters. These equations can be obtained by extension 
of the results for a single helicopter with single-cable 
suspension. Results are also given for the special case 
in which all attachments are at the e.g.’s of the rigid 
bodies; these include the equations for the point-mass 
approximation of the system. 

System and constraints- The n-body system 
shown in figure 7 consists of a load D n sus- 
pended , pendant-style, from m = n — 1 helicopters, 
SI, ... , 5m, by cables C 1 , . . . , Cm attached to the he- 
licopters at Rl, . . . ,Rm. The load and short cables 
and the hardware attaching the load to the pendant 
vertex at Rn are considered a single rigid body. If the 
cables are inelastic, then each imposes one holonomic 
constraint on the motion of a helicopter relative to the 
load by fixing the distance 

tj{r) = |Rjn| 

= \Rn* N +T N , n Rn*n n -Rf N -T Nj Rj'jj\ = toj 

(142) 

For inelastic cables 

c = 77i — 1 and d = 6n — m (143) 

Generalized velocity coordinates The general- 
ized velocity coordinates u, are to be selected so as to 
separate motion caused by cable stretching from mo- 
tion with invariant cable lengths. This is readily done 
by extension of the treatment of a single-cable suspen- 
sion from a single helicopter given in appendix B. The 
coordinates u consist of the c.g. velocity of a refer- 
ence body, taken as the load in this case, Vn * N , and 
the rn cable velocities given by their components in 


cable axes, {V jn c j j = These are aug- 

mented by the angular velocities of all n rigid bodies to 
obtain 6n coordinates for the elastic suspension. For 
the inelastic suspension, the cable-stretching coordi- 
nate £j of each cable velocity is nulled, leaving 6n — m 
coordinates. 

To assemble the kinematic relations u(u),u(u) note 
that the c.g. position of each helicopter is given by 

Rj* N = Rn* N +Rn*n N -Rjn N -Rj*ju j = 1, . . . ,m 

(144) 

whence its c.g. velocity is given in terms of coordinates 
from u as 

Vj* N = Vn* N - T NiU S(Rn*n n ) un n - T NyC j V jn cj 

+T Ntj S(Rj r jj) ujj j = 1, . . . , m (145) 

This can be rearranged to give the cable velocities in 
terms of the coordinates of v as 

Vjn cj = {Vn* N - Vf N + T NJ S(Rfjj) uj 3 
-T Nfn S(Rn*n n ) um n ) j = 

(146) 

The kinematics u(u),u(u) can be obtained from equa- 
tions (145) and (146), and are presented in part a 
of table 16. The coefficient matrices A,A~ l contain 
3(rz 4- 1) rows from the unit matrix, corresponding to 
coordinates that are in both and the remaining 
rows contain only transformations and cross products 
from Coriolis velocities. 

External forces and inertia coupling- The ex- 
ternal forces and moments on the configuration, /o, 
are assembled in part b of table 16. Its vector el- 
ements are denoted (F01, . . . ,M0n) for convenience 
in later equations. The fo term contains forces 
and moments on the rigid bodies due to weight, 
(mlg, g), aerodynamics and helicopter rotor 

output, (FA1, . . . ,MAn), and the inertia coupling 
terms X + DAu. The submatrices of A require only 
time-derivatives of transformations, and the corre- 
sponding terms in fo due to Au are all Coriolis accel- 
erations of the cables ejej x Vjn, and centrifugal ac- 
celerations of the c.g. attachment-point moment arms 
for all the bodies, a ;j x (u>j x Rj*j). 
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Suspension forces- The resultant forces and mo- 
ments applied by the suspension to the rigid bod- 
ies of the configuration fc , is given in part c of ta- 
ble 16. Since the number of cables and constraints is 
equal, the configuration vectors for the cable tensions 
{H 1, . . . , Hm} comprise a basis of the suspension force 
space which can also be used to give fc for inelastic ca- 
bles; as a result, only the cable- tension formulas differ 
for the elastic and inelastic cable models. For inelas- 
tic cables, it is necessary to invert the m x m matrix 
H t D~ l H . 

Simulation equations- Finally, the equations for 
u are assembled in part d of table 16 from u = A” 1 sf. 
Here, sf refers to the specific forces and moments 
on the configuration due to fo -f fc. Its vector el- 
ements are denoted SF1, . . . , SMn for convenience. 
The equations for Vn* N and {uijj , j = l,...,n} are 
identical to those obtained from the Newton-Euler 
equations; the cable- velocity equations depend on dif- 
ferences of specific forces and moments. 

Special case of c.g. attachments on point- mass 
bodies- A simpler set of equations is obtained if all 
cable attachments are assumed to be made at the he- 
licopter and load c.g.’s (Rj*j = 0, j = l,...,n). 

Results are listed in table 17 and are obtained by im- 
posing the c.g. attachment-point condition on the re- 
sults given in table 16. The matrices A, A -1 become 
block diagonal for the linear and angular velocity co- 
ordinates; this accounts in part for the separation of 


these coordinates into independent subsystems. Fur- 
thermore, since the suspension can apply no moments 
to any body (MCj = 0, j = the attitude 

dynamics of each rigid body are independent of all 
other coordinates of the configuration motion, except 
through the aerodynamic moments, and are defined by 
the usual Euler equation for independent rigid bodies 
(table 17). 

The remaining equations in table 17 govern the lin- 
ear velocity coordinates and are also simulation equa- 
tions for the point-mass approximation of the system. 
The c.g. velocities now depend only on cable velocities 
in equations (145) and (146), and in the assembled ve- 
locity kinematics given in table 17. The applied forces 
have no coupling with the rigid-body angular veloci- 
ties through Au, although such a coupling can occur 
in the aerodynamics. The interaction forces for the in- 
elastic suspension are now independent of the applied 
moments on any body. The corresponding cable ten- 
sions still require the inverse of an m x m matrix. This 
is given in table 17 as a sum of a diagonal matrix of 
mass ratios and the positive semidefinite Grammian 
of the cable directions, and can be inverted without 
computational difficulty. The final result is almost as 
simple as the elastic suspension formulation using rigid 
body coordinates. It uses cable-velocity coordinates 
and accommodates both inelastic and elastic suspen- 
sions with explicit calculation of cable tensions in both 
cases. 
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Table 16. Simulation equation summary: multilift with pendant suspension 


(a) Configuration kinematics: 


v = A u 



Aj,n-rj — ^NJ 3j) 3 — 1, . . . , 771 

^4i s 2n = -T Ni n S(Rn*n n ) 

= ^cj^N Aj y n+j j — 1 , . . - , 771 
Bj,2n — '^cjjN ^ 1,2 n j = 1 , ...» 771 
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Table 16. Continued. 


(b) Applied forces and inertia coupling: 

/ F0l n \ / ml g N -|- FAIn - ml A\u \ 


fo — f9 + fa — X — DAu 


FOmx 


FOrtw 


M Oli 


mm (7 ;v + FAmTv — mm A rn u 


mn gs 4- FAnTv 


MAl\ — S(wli) Jl Ldl\ 


\ M0n n J \ MAn n — S(un n ) Jn ujn n / 


where: 

Aj = j^row of A 

Aj u = -T N>cj S(wcj C j) V jn c j - T nj S 2 (ujj) Rj*jj 4 T NtU S 2 (um n ) Rn*n n j = 1, m 
MCJcj (tPcji @cj COS 4>cji 0 CJ sin 4> c j) j = l,...,m 
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Table 16. Continued. 


(c) Suspension forces: 

! fccljv 0 
0 fcc2jV 


fc = HTC 


0 


0 


-feeler — fcc2jv 
eiii o 
0 <£22 2 


0 \ 
0 


kcrriN 
— fccmjv 
0 
0 


/ TCI ^ 
TC2 

KtCtjiJ 


0 0 ... (mm n 

V -£nl„ -£n2„ -£nm n / 

O'ij = (Rj‘j x kcj)j = S(Rj*jj) Tj,N kcj N 
£nj n = (Rn*n x kcj)„ = S{Rn*n n ) T„ :N kcj N 

' (max{0, Kj {tj - toj)}, j = l,„.,m) T (Elastic cables) 
k -[i/ T D -1 if] -1 // T D~ l fo (Inelastic cables) 


) j = 1 ) •••! m 


TC = 4 
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Table 16. Concluded. 


(d) Simulation equations: 



/ 5Fljv \ 


^ (FOliv + FCl N )/ml \ 

sf = D- 1 (fo + fc) = 

SFun 

SMh 

= 

(FOnjv + FCriN)/mn 
Jl~ l {M0h + MCU) 


V SMn n / 


\jn~ 1 (M0n n + MCn n )J 


u — A~ l sf 

Vjricj — T c j t ^[SFns ~ SFjw -I- Aj,n+j SMjj + Ai^n 5Mn n ] j = 1, ...,m 

Vn*^ = SFtin 

ujj = SMjj j = 1 , n 
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Table IT. Simulation equation summary for c.g. attachments: 
multilift with pendant suspension 

c.g. Attachments: Rj*j = 0 j = 1 

Angular velocity coordinates: 

vjj = Jj l WAjj - S(wjj) Jj ujj] j = 1 , n 


Linear velocity coordinates 
(a) Configuration kinematics: 



(b) Applied forces and inertia coupling: 
/ F01^ \ 

fo = 


V FOnw / 

FOJn = rnj g N + FAj N — mj T NyCj S(u>cj C j) Vjn cj 

ucjcj = (4>cj> 0 C j cos (f> cj 1 —0 C j sin <j) cj ) 

FOun = mn + FAn n 


II 
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Table 17. Concluded. 


(c) Suspension forces: 





( FCIn \ 


/ k&N 

0 ... 

0 \ 


FC2 n 


0 

/cc2at - - • 

0 

/c = 

FCmjv 

= HTC = 

0 

0 ... 

kempt 


V TCu n / 


V — &cl/v 

"~kc2 jY i « . 

—kcm jsT j 


/ TCI \ 
TC 2 

\rcmJ 


FCj N = fccjiv FCj 
FCriN = TCj 


TC= < 


' (max{ 0, Kj (£j — toj )}, j = 1, m) T (Elastic cables) 

< 

k — [F r D -1 F] -1 F r D _1 /o (Inelastic cables) 

j kclft (FOljsr/ml - FOn^ /mn) ^ 

H t D~ l fo = 

\kcmjj (F0m^/mm - FOn^/ynn) J 


H t D- 1 H = ^[diag{^, j = 1, m] + K T K] 


K = [/ccljv, kc 2 w, ... , fccm^v] 

(d) Simulation equations: 

/ FFljv \ ^ (FOliv + FCl/\r)/ml ^ 

s/ = D- 1 (/o + /c) = 

V S Fn tv / 


\ (FOtin + FCn^)/mnJ 


u — A 1 sf 
Vn* N = SFnjsj 

Vjn cj = T CJ - at (SFn;v - 5Fjv) j = 1, •••, m 
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APPENDIX F 

LINEARIZED EQUATIONS FOR DUAL-LIFT 
SYSTEM WITH SPREADER BAR 


INTRODUCTION 

Linearized equations of motion (LEOMs) are given 
for the dual-lift system with spreader bar for static 
equilibrium reference flight conditions and assuming 
inelastic suspension. These are derived from the non- 
linear EOMs, equation (14); LEOMs for general slung- 
load systems and for general reference flight conditions 
are given in table 4, and these are specialized to static 
equilibrium in equation (41). The dual-lift system ge- 
ometry and details of the nomenclature arc given in 
appendix D (fig. 6), and notation for the configuration 
vectors and matrices in the LEOMs follows that estab- 
lished in section 6. Standard linearized aerodynamics 
are assumed for each body axis as stated in section 6. 
Various secondary effects are neglected, including in- 
terbody interferences. 

Coordinates and kinematics- Coordinates are 
often selected for linear analysis to exploit any nat- 
ural decomposition of the perturbations into nearly 
decoupled subsystems of forces, controls, and motion 
variables. For a single-symmetric aeronautical body 
aligned with the air-velocity vector, motions in the 
plane of symmetry containing (Va, ib, kb) are nearly 
decoupled from its lateral-directional motions. For the 
multibody dual-lift system in hover, all mass lies in or 
near a plane defined by the spreader bar and local ver- 
tical (i3, kN), and the system modes of motion in this 
plane are nearly decoupled from its yawing and pen- 
dulum modes of motion lateral to this plane. These 
two decompositions can be made coincident at hover 
by aligning the helicopters in or perpendicular to the 
plane of the suspension; this is done in references 6 
and 7 where a comprehensive description of the nat- 
ural modes is given. In static equilibrium at cruise 
speeds, the spreader-bar heading may be at an angle 
to the flight path, with the helicopters arranged as in 
formation flying, and with the load-bar triangle swept 
back to moderate angles from the local vertical, if there 
is significant load specific drag (fig. 8). In this case, 
the aerodynamics can couple the motion in and lateral 
to the plane of the load and spreader bar. 


For the present work, the coordinates are selected 
to maximize the longitudinal-lateral aerodynamic de- 
coupling; this is expected to be a satisfactory starting 
set for analysis of the natural modes. The generalized 
velocity coordinates u, and the kinematics v(q, u) are 
given in part a of table 18. These differ from those 
in appendix D in the use of body axis coordinates of 
the reference-point velocity relative to the mean wind, 
VA1*, and tether cable velocities relative to body axes, 
AV13 c i, AV24 c2 , where 

AV13 = V13 — Lol x R13 

AV24 = V24 - uj2 x R24 

These cable velocities and the cable axis transforma- 
tions are given in terms of cable angles relative to he- 
licopter body axes as given in part a of table 18, and 
discussed in appendix A. Note that u contains all 24 
velocity coordinates of the system with elastic suspen- 
sion. This allows an easier statement of results. Re- 
sults for the 20 coordinates of the system with inelastic 
suspension are obtained by deleting columns 6, 9, 10, 
and 12 from A corresponding to the scalar coordinates 
AV13 • kcl, AV24 • kc2, VT • i3, VT • kt of u t or 
by deleting appropriate rows and columns from the co- 
efficient matrices of the LEOMs derived from A for the 
system of 24 coordinates, as is done at the end of this 
analysis. 

The generalized position coordinates q and the kine- 
matic relation u(q,q) are included in part a of ta- 
ble 18, for the 20 coordinates of the inelastic system. 
The coordinates u, q selected here are reference-point 
coordinates consisting of the c.g. position, and ve- 
locity of body Bl, fVAl*dt , VA\\, plus additional 
position coordinates q = {Aacl, . . . , o4}, which are 
fixed in static equilibrium, and velocity coordinates 
u — {AV13 • icl, . . . ,0^44}, which are zero in static 
equilibrium. The cable angles Aacl,Aac2 in q are 
angles relative to helicopter body axes. The reference- 
point position coordinates, jVAl\dt in <7, have no 
useful physical meaning, and they do not occur in the 
LEOMs. Their derivatives, VAl*, do occur, and these 
are meaningful and convenient in aeronautical work. 
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Linearized EOMs for static equilibrium The 
LEOMs are listed in part b of table 18, for both the 
state equation for ( Sq , Su) and the second-order differ- 
ential equation for Sq. These repeat equations (46b) 
and (47), which were obtained by specializing the 
equations in table 4 for general slung-load systems 
and general (accelerating) reference flight conditions to 
static equilibrium and for reference-point coordinates, 
q , u. The primitive configuration vectors, matrices, 
and perturbations va , fba , T, Fb$, Fb va , Sq, Su, Su , 
A are defined in section 6. The stability and control 
derivatives Fb&,Fb va are assumed to be given, along 
with the reference flight condition and the correspond- 
ing values of va a , fba 0 ,T 0 . The matrix, is defined 
in table 4 and has one nonzero column for static equi- 
librium and reference-point coordinates (eq. (46a)). 

Results for the coefficient matrices of the state equa- 
tion, M 0 , Q&, Q r u , Q q , for the dual-lift system are given 
in parts c-h of table 18. These matrices are expanded 
to a working form; that is, they are partitioned into 
submatrices corresponding to the subdivision of the co- 
ordinates q, u in part a of table 18, and the submatrices 
are given in terms of natural vectors and matrices us- 
ing coordinate frames and formulations as they would 
be found in a working nonlinear simulation based on 
appendix D. MACSYMA was used to carry out routine 
matrix product expansions. 

Acceleration coefficient matrix- Results for the 
symmetric acceleration coefficient matrix, 
M 0 = A°DA 0 are given in part c of table 18. 
The matrix M~ l serves the general function of 
mapping perturbation forces and moments from 
(■ Q' u Su + Q q Sq -f Qt, A) to perturbation accelera- 
tions Su. Details of the map describing the sources of 
effects on each element of Su are inaccessible analyti- 
cally but can be given from numerical inverses of M 0 
in specific cases. 

Control term- The control term 
QsA = T 0 Fb$A is given in part d of table 18, 
and contains all control derivatives. Controls for the 
two helicopters 51,52 are represented. Load controls 
are omitted but can be added routinely if these are of 
interest. The matrix T 0 transforms the derivatives of 
the body axis components of force in Fb$ to derivatives 
of their inertial components, and the matrix A 0 does 
bookkeeping on the coordinate frames and takes cross 
products with moment arms to generate the perturba- 
tion moments due to perturbation forces. The vector 
elements of Q$A and the sum, Q f Su + Q q Sq 4- Qs A 
are perturbation forces and moments with units and 
coordinate frames consistent with the usual listing of 


forces and moments in the configuration force vectors 
fa, fc, etc. 

Velocity coefficient matrix- The velocity coef- 
ficient matrix Q f u in part e of table 18, contains 
all the stability derivatives of the system, along 
with terms in the reference velocity. The matrix 
F va j v = T 0 Fb va Tj transforms Fb va to gradients 
of the inertial components of the forces relative to the 
inertial components of the air velocity vectors. Its sub- 
matrices are conveniently indicated by the notation 
FAINvn, . * • i M A44 w defined in part e of table 18. 
The matrix A 0 in Q' u does bookkeeping on coordinate 
frames and applied moment arms. Every submatrix 
of Q ; u maps an input perturbation from a vector ele- 
ment of Su to an output perturbation force or moment 
corresponding to the elements of fa, and has units 
and input-output reference frames consistent with this 
function. Every column is associated with an element 
of Su, as indicated in part e of table 18, and this de- 
fines the input reference frame and input units for all 
its submatrices. Every row is associated with an el- 
ement of Qs A, and this defines the output reference 
frame and output units for all its submatrices. 

Position coefficient matrix- The position coeffi- 
cient matrix Q q is given as a sum of two parts, Q\ q 
(part f of table 18), which contains the reference tra- 
jectory forces fg, fba Q ; and Q2 q (part g of table 18) 
which contains a combination of stability derivatives 
and the reference velocity. The nonzero submatrices 
of Q q map angle perturbations from Sq to perturba- 
tion forces and moments corresponding to the vector 
elements of fa, and their units and output reference 
frames are consistent with this function. Each column 
of Q q gives the effect of an element from Sq; the refer- 
ence point position perturbation has no effect on the 
system motion, and all other elements of Sq have a 
nontrivial effect. 

The position coefficient matrix requires the gradi- 
ents of several configuration vectors containing A{q), 
T{q), and these entail gradients of coordinate trans- 
formations and the matrix BT. A general formula 
for the gradient of coordinate transformations is given 
in appendix A. A treatment of the gradient of BT 
is appended in part h of table 18. This matrix oc- 
curs solely in the submatrix A47 of A and contributes 
to the seventh row of Q l q in the derivatives of //I3 , 
where ju I3 is a moment due to load external forces, 
F4 = m4 g + FA4, and is given in F3 coordinates by 

/zl 3 = BT t T t N F4a/- 

= (IK J 3 ' + kN k3' sec0 3 ]R 3 *5 x F4) 3 
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A convenient scalar form is given in part h of ta- 
ble 18. The unit vectors {j3'.,k3'} are associated with 
the spreader bar, and their coordinates in JF 3j F N are 
noted in part h of table 18. The vector, //I 3 is seen 
to be a function of the variables {0*, 03 , 03 , 0a} from 
q , and the required derivatives are obtained routinely, 
and are given in part h of table 18. 

For static equilibrium, the reference trajectory forces 
in Ql q , satisfy 

FI + F2 + F3 -f F4 = 0 

where Fi is the total applied external force on Bi 
(eq. (38)). This has been used to simplify several sub- 
matrices of Ql q - The load and spreader-bar aerody- 
namics occur only in Qlu.Qlxs in columns 7 and 8, 
corresponding to <5 q 3, <5a4. In sufficiently low-speed 
flight these can be neglected, and most terms in 
columns 7 and 8 drop out, leaving only load-bar weight 
effects. An approximate view of the equilibrium forces 
on each body is given by assuming that the load- 
bar aerodynamic forces are dominated by drag, and 
that each helicopter supports its share of the load-bar 
weight and aerodynamic force: 

F3 = m3 g -I- FA3 ~ m3 g - D3 iva 

F4 =: m4 g 4- FA4 ss m4 g - DA iva 

L = F3 + F4 « (??7 3 + m4)g — (D3 4- DA) iva 

FI = ml g + FA1 = — FCl - - TCI kcl « -pL 

F2 = m2 g + FA2 = ~FC2 = -TC2 kc2 » — (1— p)L 

where iva is the direction of the reference air- velocity 
vector, VAO; L is the load to be supported by the two 
helicopters, and is transmitted to them by means of 
the tether cable forces, FC1,FC2; and p is the load 
fraction supported by the helicopter Bl. If these are 
substituted in Q 1 Q , then only the load-bar weight and 
drag are present. In general, cable forces FC1,FC2 
can include mutually cancelling tugging of the heli- 
copters against each other through the spreader bar, 
but is not present in the optimized reference config- 
uration assumed in the expressions above. The net 
external forces on the helicopters are directed oppo- 
site the cables In equilibrium. Finally, it is noted that 
all moments MAI, . . . ,MA4 drop out of the gradi- 
ent in Ql q . The terms in Ql q 6q are perturbation 
forces and moments. A typical perturbation force is 
illustrated by the contribution from the submatrix in 
row 2, column 5: 

T c i,nS(F1n) VFljv Sal — (FI x dal) c i 
dal = il <50i + jl'601 + kN <50i 


FI ~ — pTj 

This is seen to be the cable- axis components of the 
cross product of the net external force on Bl with 
the total attitude perturbation of Bl. The load-bar 
drag at cruise speeds appears in L. A typical pertur- 
bation moment is illustrated by the contribution from 
the submatrix in row 5, column 5. 

-A\ h S(F1 n ) W\ n 8a 1 = — (Rl*3 x (FI x dal))! 

which is the moment about the c.g. of Bl of the 
perturbation force described above applied at the 
tether cable attachment point, R3, and given in T\ 
components. 

The second term, Q2q , contains the stability deriva- 
tives of all bodies except the reference body, B 1, 
in combination with the reference air-velocity vector, 
VAO, and is zero in hover. More generally, its terms 
can be rationalized by forming the vector Q2 q 6q. In 
the result, all terms have the same form; for example, 
the third element is 

T c2)iV FA2N vn S(VA0 n )[W1 n Sal - W2 N 6a 2] 

= T c 2 ,7V FA2Nvn (VAO x (dal — da2))jv 

~ T C 2,n F A2Nvm 

x [kN N VAO <5(02 -0 X )- jh N VAO <5(0 2 - 0i)] 

where dal,da2 are attitude perturbations of Bl,B2 
as above, jh is perpendicular to the reference veloc- 
ity in the horizontal plane, and the expanded ap- 
proximate result assumes that the body and level- 
heading axes T\,T 2 ,Th are nearly parallel. The 
result is an aerodynamic perturbation force due to 
the perturbation velocities kN V L40 <5 (0 2 — 0i), 
jh VAO 5(02 — 0i )♦ All terms in <32^ Sq similarly con- 
sist of aerodynamic perturbation forces and moments 
from bodies B*, i = 2,3,4 due to the perturbation 
velocities kN VAO <5(0* — 0i), jh VAO <5(0* — 0i), 
which are perpendicular to the flight path and depend 
on differences in attitude perturbations from those of 
the reference body. 

Second-order ODE- Finally, the state equa- 
tion is reduced to the 20 coordinates of the sys- 
tem with inelastic suspensions by deleting the rows 
of M 0 ,Q r u ,Q qi Qs and the columns of M 0 ,Q f u cor- 
responding to the cable-stretching coordinates of u 
(scalar coordinates 6, 9, 10, and 12). The second-order 
ODE for Sq is obtained from this result as defined in 
part b of table 18; that is, M 0: Q' U are post-multiplied 
by the block-diagonal matrix, U (part a of table 18). 


110 



This modifies the columns of these matrices to convert of <5u, 6u to the Euler-angle rates of bq, bq . 
from the linear cable velocities and angular velocities 
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Table 18. Linearized equations of motion for dual-lift system 


(a) Coordinates and kinematics: 




( TV,1 

0 

0 

0 

0 

0 

0 

0 

/ FA1*! \ 

vt n 


TV.i 

TV, cl - 

-TV, C 2 

0 

^25 

^26 

^27 

0 


AV’13 cl 

* ^ 
CO 


TV,1 

TV, ci 

0 

0 

^25 

0 

-4-37 

0 


AF24 c2 

V4* n 

= A u + w = 

?V,1 

TV, ci 

0 

TV,* 

^25 

0 

A47 

-448 


VT t 

ujIi 






J 





Wli 

u;2s 







/ 




w2 2 

1^33 








I 



ui3 3 

V u 4 4 / 


V 







I / 


K w4 4 / 


Cable velocity coordinates: 

T'ci.n — Ei(Acf) c i) E 2 {A9 c1 ) Ti t N 
T C 2 ,n = Ex{A 4> c2 ) E 2 (Ae c2 ) t 2>n 


Submatrices of A: 

a 25 = -t Ni i s(m*ii + n fccii) 

A 26 = Tat , 2 S(R2*2 2 + i 2 kc2 2 ) 


™ci ^ci.jv ^13^ — AF13 c i — (R13 x wl) c i A 27 = —Tn,z 5(i?343) 

K24 c2 = T c2 , n R24 n = AV24 c2 - (R24 x w 2) c2 A 37 = T Nt3 S(R3*3 3 ) 

AV13 cl = T clil f?13i = (£1 AO c i cos A^ c i, -tl A4> cl jl) T A 47 = A 37 + T Ntt BT 
AV24 cl = T c2<n R24 2 = (t 2 A8 c2 cos A <p c2 , -t 2 A<p c2 , h) T A 48 = T NA S(R4*5 4 ) 
VT t = (fax, -£3z <j> t , (3 z) t 
BT = -S(R3*5 t ) [j3’ t j3'l + kN t k3' 3 T sec0 3 ] 


Notes: Scalar elements of BT are given in appendix D, table 14. 


WO is the mean wind. 



Table 18. Continued. 


Generalized position coordinates (inelastic suspension): 




/VA1\\ 

AV13 • icl 


AV13 • jcl 


Adel 

A V24 • ic2 


Adc2 

AV24 • jc2 


4 >t 

VT.jt 

= U(q) q = 17( 9 ) 

dl 

col 1 


d2 

co2>2 


d3 

co3s 

1 


\ d4 / 

C0A4 / 




q = U~ l u 

Aacl = (A0 c i, A8 c i) T 
Aac2 = (A$ c2 , A0 c2 ) t 

U(q) = diag{J, U 22 , U 33 , -£3z, TTli, 1T2 2 , WS 3 , W4 4 } 

U~ 1 (q) = diag{I, U 2 ~ 2 \ U 3 ~ 3 l , -1/tZz, Wl^\ W2,\ W3 3 \ W4, 1 } 


U 22 = 


0 £1 cos A0 c i 

— £1 0 


££33 = 



£2 cosA^^' 
0 


Note: {Wii, Wir 1 1 = 1,2, 3, 4} 


defined in appendix A, equations (67) and (68) 
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Table 18 . Continued. 


(b) LEOMs for static equilibrium: 

8q = U~ l 8u 

$ii = M 0 1 (Q' u 8u + Q q 8q + Q& A) 

or 

M*o Sq -Q'q 8q- Q q 8q = Q & A 
where 

M 0 = A T 0 D A 0 

Q'u = A T a F vaN A a - Ag D & 

Qq = Qlg + Q2q 

Ql g = [VJ A T (q) (j fg + fa 0 ) + Aj V T q T(q ) /6a 0 )] g „ 
Q2, = A T a F vaN [V^ (T 0 T(qf (v a -w)+ A(q) «„)],„ 
Q* = A T 0 T 0 Fb s 
FvaN = To Fb va TJ 

£ u : defined in equation (46) for static equilibrium, 
and 

M ' 0 = M 0 U 0 
Q'q = Q'u Uo 
Uo = U(q a ) 
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Table 18. (continued) 
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Table 18 . Continued. 


(d) Control vector: Q s A = A? T 0 Fb 6 A? 





/ T ltN \ 

0 



T c l,N 

0 



—Tc2,N 

0 



0 


Ai + 



VJ X MAU 



^25 

0 



^26 

0 



^7 

Vo/ 


V 

V 0 y 


Tn, 2 Vj 2 FA2 2 


Note: Ai, A 2 — control perturbations for helicopter #1,#2. 


(e) Velocity coefficient matrix: Q' u — A J F va H A 0 — Aj D 

6VV, 6VU cl <5^24^ SVT t 6u>\ } 6a;2 2 5o;3 3 
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t- 
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\ 
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2: 

e5 

* 

«*1 

> 
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hM 

<o 

o 

c! 

+ 

t- 
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o 
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<o 

Cj 

CO 

o 

2: - 

o 

O 

50 

CO 

o 

o 

5: 


o 
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o 

L s 

et 

e5 

K 

i 





^ 8 

i o' 2 O t-rl 
^ ^ >!*• ^ 
c* C o < to I 

< <N 3 •* a 

i- E 2 S - 

*+~ + 3 

E? 5 SI 


o o 

* h.% 

E7 ^ 


N V M ' C* 7 

^ ^ ^ + 

Jo _.£ x> 

o o o 


s $ 


2: ^ ^ ^ ^ 

cT c5 > £ r* + 

fc. fa. £ > S - 3 

— •*■ CM > . > 


2: ^ * 
>" « o CM CM TT 

tx, &. *&, 
.> £ ^ + 


< > s j j 3 j 
g ^ 4- ^ ^ 5 ^ 

C*. ^ 

II II II II II + If 

CO M* O CO b- X 

fioooo so 


o 

y 

f? 

5:' 





t? 

6? 



E? 




~cr 

CO 

CM 

Cx, 

TT 

CO 

CM 

£ 

CM 


-cf 

CO 

e 

E? 

’o 

E? 

4J 

E? 

Q 

Cs 


Ct 

Q 

CJ 

cS 

c! 

2 

E7 

2: 

e5 

2: 

ci 

eS 

**4 







r-4 

- 

CO ^ 


r 2:* 



E-i 

£ S 


<3 

^ rr 




^ § " II II II II II 

£ Cl P ~ o ~ f- - 

**-< rn tn r- h- ct 

Q Q S CoO OOO 
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Css = MAMu + MAAAvn Am + A{* Cig c '* = £>^234 A 2 s + Jn 23 i SVO 

Note: j 4 2 7, ...» A 48 , m 234 , iW “]7 defined in parts a and c of table 18 



(f) Position coefficient matrix: Q q = Ql,, + Q2 q \ Q\ q = \V q (A(q) T ( fg + fn 0 ) + A„ T(q) fba 0 )] q0 : 

6 [ VAX* i foci foc2 fol 6a2 8a3 5a4 



£ 

CO 



5; * £ 

N M 

^ ^ 


(t-j ^ 

t ^ ^ - 

^ ^ ^ 

CO Co Co > 

!l il II II 



* 

* 




cq 

It 


n 




o> 


118 

1 III 


^ 25 i * « • > ^48 defined in part a of table 18; /AI 3 and its derivatives given in part h. 
Scalar elements of Wb^ } iCluiC2 2 given in appendix A, equations (76) and (94). 



Table 18 . Continued. 
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Q2 77 = [MA33 vn + Aj 7 FA3 N vn ]S(VA(\ n ) W3 n 
Q 2sa = [MA44vn + F A4Nv n\S(VA0^/) W4 jv 

Notes: FA2N V n, • • • ,MA44 VN defined in table FI (<1). /W--,Ais defined in part a of table 18. Wb^, 6 = 2,3, 4 given in appendix A. 



Table 18. Concluded. 


(h) Derivatives of /il 3 : dfil^/d<f> ti Vj 3 fi 1 3 : 

jj,l 3 = BT T T t)N F4 n = -[j3 3 j3'^ 4 A:3 3 kNjj- sec0 3 ] S(F4 N ) T Nit R3*5 t 

where: 

j3' = — sin iN 4 cos^ 3 jN = cos <fi 3 j3 — sin 0 3 k3 
k3' = sin 4> 3 j3 4- cos 0 3 k3 

then: 

dfjl 3 /d<pt — Zl 3 

9/il 3 /003 = [ 3^3 kN% sec 0 3 + k3' 3 j3' N T ] jj2 N 
dfil 3 /d0 3 = Z2$ — sec 9 3 tan(9 3 A;3 3 kNj^ ji2^ 
duU/dth = Zh + jZ'i K" t n p2 N 

where: 

i3" = cos ^ 3 iN 4 sin?/? 3 jN 
(Fix — S(F4^) T/v,* R3*h t 
BTl = j3' 3 j3' n t + k3' ?j kNj, sec 9 3 
[Zl 3 , Z2 3 , Z3 s ] = BTl S(F4jv) T N>t S(R3*5 t ) Wt t 
and: Wt t = [it t , j3' t , kN t ] (appendix A, eq. (77)) 
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